προβλημα σε εξίσωση του matlab

[theo2008]

Καλησπερα σε ολους.Εχω ενα προβλημα με μια εξισωση στο matlab και χρειαζομαι την βοηθεια σας(οχι την επιλυση φυσικα).Ειχα στειλει και αλλες φορες την ιδια απορια απλα τωρα στελνω ολοκληρωτικα αυτο που θελω να ρωτησω και που εχω κολλησει.

Εχω την εξης συναρτηση.

Cdot=qAL∂/∂v∫(1/(1+(exp(E-qVdot)/kT))*(gNdot/sqrt(piM/2))*exp(-2((E+Edot+qVdot)/M)^2)dE.

 

Αυτη η εξισωση αναφερεται σε θεωρητικα μοντελα αναλογικων στοιχειων και δεν επιλυεται αναλυτικα.Πρεπει να υλοποιηθει αριθμητικα.

 

Αυτην την εξισωση πρεπει να την περασω στο matlab.Οι παραμετροι που εχω ειναι:

A=12.56*10^-4;,L=0.44;,G=2;,T=300;,q=1.6*10^-19;,k=8.62*10^-5;,

Edot=0.01;,Ndot=10^8;,M=0.001;

Το Edot μπορει να παρει τιμες 0.01-1

το Ndot 10^8-10^15

το M 0.001-0.5

Το ολοκληρωμα θελω να ειναι απο (0,inf)

Το V εχει τιμες V=-12:0.05:2;

Επισης ο τυπος για το vdot ειναι Vdot=V/L;

Με αυτην την συναρτηση θελω για καθε τιμη του V να βγάλω και μια τιμή για το Cdot.Μετα αφου βγουν αυτες οι τιμες να κανω μια γραφικη V,Cdot.

 

Ο κωδικας που εχω σκεφτει ειναι:

Πρωτα κανω ενα M-file με το ονομα integral3:

 

function total=integral3(vdot)

A=12.56*10^-4;

L=0.44;

G=2;

T=300;

q=1.6*10^-19;

k=8.62*10^-5;

Edot=0.01;

Ndot=10^8;

M=0.001;

de=0.01;

E=0:0.01:1000;

expe=((E+Edot+(q*vdot)/L)/M).^2;

den=1+exp((E-(q*vdot)/L)/(k*T));

static=(G*Ndot)/(sqrt(pi*M/2));

sinar2=exp(-2*expe);

sinar1=1./den;

sinar3=sinar1.*static;

final=sinar3.*sinar2;

subtotal=sum(final)*de;

total=q*A*L*subtotal;

 

 

 

Μετα αλλο ενα M-file:

cdot1(1)=integral3(-12);

cdot2(281)=integral3(2.05);

for i=2:281

vdot(i)=-12+(i-1)*0.05;

cdot1(i)=integral3(vdot(i));

cdot2(i-1)=cdot1(i);

cdot(i)=1e11*(cdot2(i)-cdot1(i))/0.05;

end

plot(vdot,cdot,'.:')

 

 

 

Το προβλημα μου ειναι πως σε καθε επαναληψη του κωδικα μου βγαινει το ιδιο αποτελεσμα για την Cdot με αποτελεσμα η γραφικη να βγαινει μια ευθεια που προφανως ειναι λαθος.Μπορειτε να μου πειτε που κανω λαθος ή να μου προτείνετε καποιον αλλο τροπο?Επισης αυτη η εξισωση μηπως θα μπορουσε να επιλυθει με αριθμητικη ολοκληρωση Gauss?Οχι οτι γνωριζω πολλα πραγματα για αυτην την μεθοδο..

Ευχαριστω

[theo2008]

επανερχομαι με το ιδιο θεμα που ειχα προβλημα και παλιοτερα αυτην την φορα με νεα απορια!!!Η συναρτηση ειναι αυτην που σας δινω παραπανω ή αλλιως την στελνω σε αρχειο word.Οπως αναφερω αυτην την εξισωση θελω να την επιλυσω αριθμητικα για v=-12:0.05:2.Υστερα απο προσπαθειες εχω καταληξει στο εξης.

Εχω ενα αρχειο m-file με ονομα integral που ειναι το εξης:

>
function subtotal=integral(v)
A=12.56e-4;
L=0.44;
G=2;
T=300;
q=1.6e-19;
k=8.62e-5;
Edot=0.05;
Ndot=5e9;
DE=0.25;
syms E;
denom=1+exp((E-(v/L))/(k*T));
function1=1./denom;
static=(G*Ndot)/(sqrt((pi*DE)/2));
ek8eths=((E+Edot+(v/L))/DE).^2;
function2=exp(-2*ek8eths);
final=function1.*function2.*static;
inline_final = vectorize (inline (char (final)));
subtotal = quad (inline_final, 0, 1e12);

 

και υπολογιζει το ολοκληρωμα,και ενα αρχειο m-file με ονομα derivative που υπολογιζει την παραγωγο:

>
cdot1(1)=integral(-12);
cdot2(281)=integral(2.05);
for i=2:281
  v(i)=-12+(i-1)*0.05;
   cdot1(i)=integral(v(i));
   cdot2(i-1)=cdot1(i);
   cdot(i)=(cdot2(i)-cdot2(i-1))/0.05;
end

 

Το προβλημα ειναι πως ενω στο m-file integral οταν δινω μια τιμη στο v χειροκινιτα μου βγαζει αποτελεσμα,οταν παω να τρεξω το m-file derivative τοτε μου δινει μονο ενα αποτελεσμα που προφανως ειναι λαθος μια και 8ελω για καθε v και διαφορετικο Cdot καθε φορα...

Εχετε καμια ιδεα για το τι πρεπει να διορθωσω?

συναρτηση.doc

[insomniaK]

Να είμαι ειλικρινής, δεν άνοιξα ούτε τη Matlab ούτε το αρχείο για να δω τι κάνει το script σου, αυτό που προτείνω είναι να βάλεις

 

>for i=2:281
  v(i)=-12+(i-1)*0.05;
display(v(i))
   cdot1(i)=integral(v(i));
display(cdot1(i))
   cdot2(i-1)=cdot1(i);
display(cdot2(i))
   cdot(i)=(cdot2(i)-cdot2(i-1))/0.05;
display(cdot(i))

end

 

και να βλέπεις κατά τη διάρκεια του τρεξίματος τι υπολογίζεται.

[insomniaK]

Να είμαι ειλικρινής, δεν άνοιξα ούτε τη Matlab ούτε το αρχείο για να δω τι κάνει το script σου, αυτό που προτείνω είναι να βάλεις

 

>for i=2:281
  v(i)=-12+(i-1)*0.05;
display(v(i))
   cdot1(i)=integral(v(i));
display(cdot1(i))
   cdot2(i-1)=cdot1(i);
display(cdot2(i))
   cdot(i)=(cdot2(i)-cdot2(i-1))/0.05;
display(cdot(i))

end

 

και να βλέπεις κατά τη διάρκεια του τρεξίματος τι υπολογίζεται.

[theo2008]

Σευχαριστω πολυ insomnia K.Μια χαρη ή εσυ ή οποιος αλλος μπορει επειδη μεχει κουρασει πολυ αυτη η συναρτηση που εχω να περασω στο matlab δειτε την λιγο στο αρχειο συναρτηση.doc που εχω επισυναψει και πειτε μου τουλαχιστον αν την περναω σωστα στο matlab απο το ολοκλήρωμα και μετα γιατι αυτο ειναι που με ενδιαφερει πιο πολυ.

Εγω την εχω περασει ετσι:

>
paronomastis=1+exp((E-(q*v))/(k*T));
function1=1./paronomastis;
static=(G*Ndot)/(sqrt((pi*DE)/2));
ek8eths=((E+Edot+q*v)/DE).^2;
function2=exp(-2*ek8eths);
final=function1.*function2.*static;

Νομιζω οτι ειναι σωστος ο τροπος μου!!!Εχω στειλει τοσες πολλες ερωτησεις εδω στο insomnia γιατι ειναι μια μελετη που κανω σε αυτην την εξισωση και δυστυχως μου προκυπτουν συνεχεια αποριες που δεν μπορουν να με βοηθησουν τα tutorials που εχω!!!!

[theo2008]

Σευχαριστω πολυ insomnia K.Μια χαρη ή εσυ ή οποιος αλλος μπορει επειδη μεχει κουρασει πολυ αυτη η συναρτηση που εχω να περασω στο matlab δειτε την λιγο στο αρχειο συναρτηση.doc που εχω επισυναψει και πειτε μου τουλαχιστον αν την περναω σωστα στο matlab απο το ολοκλήρωμα και μετα γιατι αυτο ειναι που με ενδιαφερει πιο πολυ.

Εγω την εχω περασει ετσι:

>
paronomastis=1+exp((E-(q*v))/(k*T));
function1=1./paronomastis;
static=(G*Ndot)/(sqrt((pi*DE)/2));
ek8eths=((E+Edot+q*v)/DE).^2;
function2=exp(-2*ek8eths);
final=function1.*function2.*static;

Νομιζω οτι ειναι σωστος ο τροπος μου!!!Εχω στειλει τοσες πολλες ερωτησεις εδω στο insomnia γιατι ειναι μια μελετη που κανω σε αυτην την εξισωση και δυστυχως μου προκυπτουν συνεχεια αποριες που δεν μπορουν να με βοηθησουν τα tutorials που εχω!!!!