partielle DGL mit Randwerten < Matlab < Mathe-Software < Mathe < Vorhilfe
|
Aufgabe | [mm] \bruch{dT_i}{dt}=\bruch{\lambda}{c*p}*\left( \bruch{A_{i+1}-A_{i-1}}{2*(L)^2}*\bruch{T_{i+1}(t)-T_{i-1}(t)}{2*L} \right)-\bruch{2*\wurzel{\Pi}*\delta*T_i ^4(t)}{\wurzel{A_i}*c*p}
[/mm]
Durch Matlab lösen lassen. |
Moin moin,
ich habe diese wunderschöne DGL, die eigentlich mal eine partielle DGL war. Ich hatte sie räumlich diskretisiert (deshalb der Index i). Dazu habe ich natürlich noch die Randwerte. Ich weiss jetzt nicht, wie ich das in Matlab eingebe. Ich habe eigentlich zwei Varianten:
1. ich diskretisiere es auch noch mal zeitlich. Das habe ich auch bereits gemacht, aber trotz einer sehr klein gewählten Schrittweite (1/6000) funktioniert es nicht. Die Werte sollten sich eigentlich stätig erhöhen, aber die sinken teilweise wieder. Und alles ganz merkwürdig. Aber das wäre eigentlich meine liebste Variante. Ist zwar relativ ungenau, aber es soll eh hauptsächlich eine grobe Entwicklung zeigen.
2. ich gebe in Matlab die DGL ein und lasse es davon lösen. Ich weiss, dass es mit bvp4c gemacht wird. Ich habe dazu auch ein Bsp. gefunden, dass ich fast vollständig verstanden habe. Aber da habe ich noch Hoffnung. Dabei bereitet mir Kopfzerbrechen, dass ich die Indexe (i, i-1, i+1) habe. Da habe ich keine Idee, wie ich das Matlab sinnvoll beibringe.
Ansonsten kann ich mal noch was zu den Randwerten sagen. Ist jetzt aber auf die erste Variante mit T_(i,k) bezogen:
[mm]T_{0,k}[/mm] und [mm]T_{n+1,k}[/mm] sind für jeden Zeitpunkt k gleich
[mm]T_{i,0}[/mm] und [mm]T_{i,m} [/mm] sind für jedes Element bekannt.
Wäre schön, wenn mir jemand helfen kann. Egal bei welcher der Varianten. eswürde eventuell schon helfen, wenn mir jemand beim sortieren der Gedanken helfen könnte. Wenn es hilft, kann ich auch noch die partielle DGL und die vollständig diskretisierte Formel posten. Wollte nur nicht noch mehr abschrecken. Also bin dankbar für alles.
MfG Wolf-im-Schafspelz
p.s.: Dringlichkeit ist so eine Sache. Ich will es natürlich so schnell wie möglich haben, aber zur Not, kann es auch etwas dauern.
edit: Soll ich etwas umformulieren oder anders beschreiben?
Ich habe diese Frage in keinem Forum auf anderen Internetseiten gestellt.
|
|
|
|
Moin moin,
ich war natürlich zwischendurch fleißig und bin glaube ich auch etwas weiter gekommen. Aber leider funktioniert mein Programm nicht. Ich bekomme hier Fehlermeldungen, mit denen ich nicht wirklich was anfangen kann. Deshalb hoffe ich mal ganz stark, dass hier jemand vorbeikommt, der Ahnung von Matlab hat.
function Hitze=Hitze()
%SIMULATIONSPARAMETER DEFINIEREN:
n=50; %Anzahl der Zylinderelemente
dt=1; %Intervall der Zeitdiskretisierung in sec
t=20; %Anzahl der Iterationsschritte (t*dt ist also die simulierte Zeitspanne)
t_pause=0.1; %Zeit zwischen dem Anzeigen zweier Simulationsschritte im Plot (t_pause=dt entspricht Echtzeit)
%PHYSIKALISCHE KOEFFIZIENTEN UND ANFANGSWERTE DEFINIEREN:
T_Start=300; %skalare Temp. in Kelvin; wird allen Zylinderelementen bis auf das erste und letzte zugewiesen
P_Heizung=600; %Heizleistung in Watt
T_Wasser=285; %Kühlwassertemp. in Kelvin
l=0.115; %Länge des Stabs in Meter
lambda=401; %80.2; %alles in SI-Einheiten
c=385;
rho=8920;
sigma=5.6704e-8;
alpha=1600; %mittlerer Wärmeübergangskoeff. Stab--Wasser
%GEOMETRIE DEFINIEREN:
d_aussen=0.02;
d_innen=0.008;
x1=0.2; %Anfang der Querschnittsänderung auf l bezogen
x2=0.6; %Ende der Querschnittsänderung auf l bezogen
xH=0.06;
xK=0.1;
n1=round(n.*x1); %Zylinderelement bei dem die lin Querschnittszunahme anfängt
n2=round(n.*x2); % " " " " endet
nH=round(n.*xH); %Anzahl der Zylinderelemente am linken Rand, auf die die Heizleistung einwirkt
nK=round(n.*xK); %Anzahl der Zylinderelemente am rechten Rand, die in das Kühlwasser eintauchen
%Querschnitt einstellen
d=d_aussen*ones(1,n);
for j= 1:(n1)
d(j)=d_innen;
end
for j= (n1):(n2)
d(j)=d_innen+(j-n1)*(d_aussen-d_innen)/(n2-n1);
end
global A;
A=pi/4 * d.^2; % durch den Punkt wird es Elemntenweise berechnet
global T;
T(1,:)=ones(1,n);
T(1,:)=T_Start * T(1,:);
tspan=0:dt:t;
for j=2:t
T(j,:)=ODE45(ODE45_RHS,tspan,T(j-1,:));
%T(j,:)=ODE45(T(j-1,:),0,dt,ODE45_RHS);
end
%plotten:
x1=linspace(0,l,n);
for i=1:t
hold on;
plot(x1,T(i,:));
xpause(1e6*t_pause);
end
%%%%%%%%% DGL: %%%%%%%%%%%
function T_dot=ODE45_RHS(t,T) %Die Funktion erwartet T als Zeilenvektor
global T A;
n=50; %Anzahl der Zylinderelemente
dt=1; %Intervall der Zeitdiskretisierung in sec
t=20; %Anzahl der Iterationsschritte (t*dt ist also die simulierte Zeitspanne)
t_pause=0.1; %Zeit zwischen dem Anzeigen zweier Simulationsschritte im Plot (t_pause=dt entspricht Echtzeit)
T_Start=300; %skalare Temp. in Kelvin; wird allen Zylinderelementen bis auf das erste und letzte zugewiesen
P_Heizung=600; %Heizleistung in Watt
T_Wasser=285; %Kühlwassertemp. in Kelvin
l=0.115; %Länge des Stabs in Meter
lambda=401; %80.2; %alles in SI-Einheiten
c=385;
rho=8920;
sigma=5.6704e-8;
alpha=1600; %mittlerer Wärmeübergangskoeff. Stab--Wasser
d_aussen=0.02;
d_innen=0.008;
x1=0.2; %Anfang der Querschnittsänderung auf l bezogen
x2=0.6; %Ende der Querschnittsänderung auf l bezogen
xH=0.06;
xK=0.1;
n1=round(n.*x1); %Zylinderelement bei dem die lin Querschnittszunahme anfängt
n2=round(n.*x2); % " " " " endet
nH=round(n.*xH); %Anzahl der Zylinderelemente am linken Rand, auf die die Heizleistung einwirkt
nK=round(n.*xK); %Anzahl der Zylinderelemente am rechten Rand, die in das Kühlwasser eintauchen
[m,n]=size(T);
dx=l/n;
T_dot=zeros(m,n);
[mm] T_dot(1)=lambda/(c*rho*A(1))*(A(1)*(T(2)-2*T(1)+T(2))/(dx^2))-(1/(c*rho*A(1)))*(2*sqrt(pi*A(1))*sigma*(T(1))^4 [/mm] -P_Heizung/(2*dx*nH));
for i= 2:(nH)
[mm] T_dot(i)=lambda/(c*rho*A(i))*((A(i+1)-A(i-1))*(T(i+1)-T(i-1))/(4*dx^2)+A(i)*(T(i+1)-2*T(i)+T(i-1))/(dx^2))-(1/(c*rho*A(i)))*(sqrt(pi*A(i))*sigma*(T(i))^4-P_Heizung/(2*nH*dx));
[/mm]
end
for i= (nH+1):(n-nK)
[mm] T_dot(i)=lambda/(c*rho*A(i))*((A(i+1)-A(i-1))*(T(i+1)-T(i-1))/(4*dx^2)+A(i)*(T(i+1)-2*T(i)+T(i-1))/(dx^2))-(2/(c*rho))*sqrt(pi/A(i))*sigma*(T(i))^4;
[/mm]
end
for i= (n-nK+1):(n-1)
[mm] T_dot(i)=lambda/(c*rho*A(i))*((A(i+1)-A(i-1))*(T(i+1)-T(i-1))/(4*dx^2)+A(i)*(T(i+1)-2*T(i)+T(i-1))/(dx^2))-(2/(c*rho))*(sqrt(pi/A(i)))*(sigma*(T(i))^4+alpha*(T(i)-T_Wasser));
[/mm]
end
[mm] T_dot(n)=(lambda*A(n)*(T(n-1)-T(n))/dx+(dx*2*sqrt(pi*A(n))+A(n))*(alpha*(T_Wasser-T(n))-sigma*(T(n))^4))/(c*rho*A(n)*dx);
[/mm]
Es ist mir klar, dass es jetzt nur schwer zu durchschauen ist. Ich hoffe, dass mir dennoch jemand helfen kann.
Also die Formel habe ich dank großer Hilfe so angepasst bekommen. Nur hat da die Person nicht wirklich Ahnung von Matlab. :-(
??? Error using ==> exist
Unknown command option.
Error in ==> [mm] C:\MATLAB6p5\toolbox\matlab\funfun\private\odearguments.m
[/mm]
On line 75 ==> if (exist(ode)==2) % M-file
Error in ==> [mm] C:\MATLAB6p5\toolbox\matlab\funfun\ode45.m
[/mm]
On line 155 ==> [neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, args, ...
Error in ==> [mm] C:\Dokumente [/mm] und [mm] Einstellungen\Visitec\Desktop\Neuer Ordner\Hitze.m
[/mm]
On line 57 ==> T(j,:)=ODE45(ODE45_RHS,tspan,T(j-1,:));
Mit dieser Fehlermeldung kann ich nicht wirklich was anfangen. Bis zum Aufruf der ODE45-Funktion klappt aber alles. So dass es wohl irgend etwas damit zu tun haben muss.
MfG Wolf-im-Schafspelz
p.s. brauche wirklich Hilfe!!! Danke schön.
|
|
|
|
|
Hallo ist da jemand der mir helfen kann?
Habe hier im Forum leider noch keine einzige nützliche Antwort bekommen. :-(
Aber zumindest konnte ich schon 2-3 Leuten helfen. aber dafür bin ich eigentlich im Augenblick nicht hier unterwegs.
MfG Wolf-im-Schafspelz
|
|
|
|
|
Status: |
(Mitteilung) Reaktion unnötig | Datum: | 17:20 Sa 29.07.2006 | Autor: | matux |
$MATUXTEXT(ueberfaellige_frage)
|
|
|
|
|
Moin moin,
naja wahrscheinlicher ist, dass sich hier niemand herum treibt, der auch mal bei Matlab nachschaut ob er/sie helfen kann.
Denn an der Zeit kann es nicht so wirklich gelegen haben. 2 Wochen sind doch eine ganz schön lange Zeit.
Zum Glück habe ich am Freitag dann doch noch in einem anderen Forum gepostet und da hat sich schnell jemand gefunden der mir helfen konnte. Ich weiss, damit habe ich gegen eure Regeln verstossen, aber was solls.
Lebt wohl und ich hoffe für euch das sich hier was ändert.
MfG Wolf-im-Schafspelz
|
|
|
|
|
Hallo Wolf-im-Schafspelz,
> Ich weiss, damit habe ich
> gegen eure Regeln verstossen, aber was solls.
Die Regeln besagen nur das man den link auf die andere Diskussion mit angeben soll. Dies kannst Du ja immer noch tun. Oder geht das etwa nicht weil es doch nicht erst letzten Freitag war? ....
viele Grüße
mathemaduenn
|
|
|
|