Rohrströmung gerades Rohr mit Plotten des Ausflussprofils

WICHTIG: Der Betrieb von goMatlab.de wird privat finanziert fortgesetzt. - Mehr Infos...

Partner:

Forum

Option

[Erweitert]

• Diese Seite per Mail weiterempfehlen

Gehe zu:

ng83

Forum-Anfänger


	Beiträge: 43

	Anmeldedatum: 26.11.12

	Wohnort: Deutschland

	Version: ---

Verfasst am: 12.01.2015, 09:02 Titel: Rohrströmung gerades Rohr mit Plotten des Ausflussprofils

Hallo,
ich möchte eine Strömung in einem geraden Rohr diskretisieren und das Ausflussprofil plotten

Meine Gleichung dazu lautet:

(u^(n+1)-u^n)/dt+(u^n grad)u^(n+1)+grad p^n=1/Re \Delta u^(n+1)

Dazu habe ich mit der Poissongleichung begonnen.
Das zugehörige fuction-file sieht wie folgt aus:

Code:

function [ u , x , y ]= poissonfd(a,c,b,d,nx,ny,f,g)
%??????????????????????????????????????????????????????????
% poissonfd loest die 2d?Poisson?Gleichung mit dem
% 5 Punkte Finite-Differenzen-Schema im Rechteck
% (a,b)X(c,d) mi t der Dirichlet-Randbedingung
% u(x,y)=g(x,y) fuer (x,y) auf dem Rand des Rechtecks .
% nx : Anzahl der inneren Gitterpunkte in x?Richtung
% ny : Anzahl der inneren Gitterpunkte in y?Richtung
% f : MATLAB?Function der rechte Seiten der Poisson?Gleichung
% g : MATLAB?Function der Randbedingung
%??????????????????????????????????????????????????????????
%%%Schrittweiten
nx=nx+1;
ny=ny+1;
hx = (b-a)/nx;
hy = (d-c)/ny;
nx1 =nx+1;
hx2 = hx^2;
hy2 = hy^2;

dim = (nx +1)*(ny +1);
K = sparse(dim,dim); %%% Matrix der FDM
rhs = zeros(dim,1); %%% Rechte Seite der FDM

x = [ a : hx : b ];
y = [ c : hy : d ];
[ X, Y ] = meshgrid(x, y);%ndgrid=lexikografische reihenfolge
rhs2 = f(X, Y)';
rhs2 = rhs2(:);
dirichletValues = g(X, Y)';
dirichletValues = dirichletValues(:);
nbound = [1 : nx1, dim-nx : dim, 1 : nx1 : dim-nx,nx1 : nx1 : dim];%Randknoten
nbound=unique(nbound);
ninternal = setdiff( 1 : dim , nbound );%innere Knoten

rhs=rhs2;

%MAtrixeintraege

e = ones(1,(nx-1)*(ny-1));
I = repmat(ninternal,1,5);
J = [ninternal,ninternal-1,ninternal+1,ninternal-nx1,ninternal+nx1];
aij = [(2/hx2+2/hy2)*e,-1/hx2*e,-1/hx2*e,-1/hy2*e,-1/hy2*e];
K=sparse(I,J,aij,dim,dim);

rhs=rhs-K(:,nbound)*dirichletValues(nbound);%Randknoten werden extrapoliert
K = K( ninternal,ninternal ); %%% Extraktion des Subsystems
rhs = rhs( ninternal ); %%% der inneren Freihei t sgrade
utemp = K\ rhs; %%% Loesen des LGS
uh = dirichletValues; %%% Setzen der Randwerte
uh ( ninternal ) = utemp; %%% Loesung an inneren Knoten

u=reshape(uh',nx+1,ny+1)';
x = [ a : hx : b ];
y = [ c : hy : d ];

Funktion ohne Link?

das Script:

Code:

x_1=0;%a
x_2=15; %b
y_1=0;%c
y_2=1;%d
f=@(x,y) 0*x+1;
g=@(x,y) 0*x;
nx=30;
ny=30;

%%% Loesung der Randwertaufgabe
[ u , x , y ]= poissonfd( x_1 , y_1 , x_2 , y_2 , nx , ny , f , g );

[ X, Y ] = meshgrid(x, y);%ndgrid=lexikografische reihenfolge
solution = g(X, Y);

figure(1)
surf ( X , Y , u);
%print ( gcf , '-djpeg100','loesungPoisson.jpg');
xlabel('x');
ylabel('y');
title('numerische Lösung');

figure(2)
surf ( X , Y , solution-u );
xlabel('x');
ylabel('y');
title('Fehler in Berechnung');
%print ( gcf , '-djpeg100','solutionPoisson.jpg');

Funktion ohne Link?

jetzt möchte ich diese zum zeitpunkt n+1 ermitteln.
Habt ihr eine Idee wie ich das realisieren könnte?

Danke für Eure Hilfe vorab.

VG ng

ng83

Themenstarter

Forum-Anfänger


	Beiträge: 43

	Anmeldedatum: 26.11.12

	Wohnort: Deutschland

	Version: ---

Verfasst am: 13.01.2015, 09:08 Titel:

Wirklich niemand, der mir sagen kann wie der n+1-Zeitpunkt aussehen kann?
Meine Idee war es eine for-Schleife zu nutzen:

Code:

for n=1:10
u(n+1)=u(n)*dt;
end

Funktion ohne Link?

und danach vllt die neuen werte ausgeben lassen.

Nur weiß ich nicht wie ich die um meinen Code basten soll.


Einstellungen und Berechtigungen
Beiträge der letzten Zeit anzeigen:	Du kannst Beiträge in dieses Forum schreiben. Du kannst auf Beiträge in diesem Forum antworten. Du kannst deine Beiträge in diesem Forum nicht bearbeiten. Du kannst deine Beiträge in diesem Forum nicht löschen. Du kannst an Umfragen in diesem Forum nicht mitmachen. Du kannst Dateien in diesem Forum posten Du kannst Dateien in diesem Forum herunterladen

Impressum | Nutzungsbedingungen | Datenschutz | FAQ |

RSS

Hosted by:

Copyright © 2007 - 2025 goMatlab.de | Dies ist keine offizielle Website der Firma The Mathworks

MATLAB, Simulink, Stateflow, Handle Graphics, Real-Time Workshop, SimBiology, SimHydraulics, SimEvents, and xPC TargetBox are registered trademarks and The MathWorks, the L-shaped membrane logo, and Embedded MATLAB are trademarks of The MathWorks, Inc.