Matlab-harjoitustyö

Seuraa 
Viestejä45973
Liittynyt3.9.2015

Tietääkö joku suomenkielistä foorumia, jossa olisi matlab-asiantuntemusta? Vai onko täällä Matlab-asiantuntijoita?

Kyse on harjoitustyöstä, jossa on ratkaistava numeerisesti yhtälö jossa tuntematon integraalien sisällä. Kysessä on 2 määrättyä integraalia, joilla on eri rajat, lisäksi ratkaistavassa lausekkeessa tuntematon jonka arvo halutaan ratkaista esiintyy kahdessa muussakin termissä, jotka eivät ole integraalin sisällä.

Kommentit (15)

Vierailija
kiltanen
Tietääkö joku suomenkielistä foorumia, jossa olisi matlab-asiantuntemusta? Vai onko täällä Matlab-asiantuntijoita?

Kotimaista matlab-tietämystä voit kysellä esim. usenet-ryhmistä sfnet.tiede.*

Mihin homma tyssää siinä integraaliyhtälön numeerisessa ratkaisemisessa?

Vierailija

Katso ratkaistavaa lauseketta osoitteesta http://www.lce.hut.fi/teaching/S-114.14 ... itus_1.pdf

kyseessä on kaava (6). Kurssin sivuilta löytyvän aputiedoston pohjalta olen rakentanut *.m-tiedoston, joka tulostaa näytölle me=9.

Paljastamatta koko koodia laitan todennäköisen ongelmakohdan tähän. EF:n ratkaisemiseen käytetään fzeroa, jonka kutsu on jäljempänä.

x=quadl(@integ_ele, EG*e, 10000,[],[],E, T, EF, EG);
y=quadl(@integ_auk, -10000, 0,[],[], E,T,EF,EG );
%******************************************
% Lauseke, jonka nollakohta etsitään
%------------------------------------------
function arvo = nollattava(E, T, EF , EG);
% Sijoita lausekkeen arvo muuttujaan arvo.
%
% Elektroni- ja aukkotiheydet sekä ionisoituneet akseptori- ja donoritiheydet
% tehtävän annon kaavasta (6) muokattuna.
r=NA*s/(4*exp((EA-EF)*e/(k*T))+1);
o=ND*s/(2*exp(-(EG-EC-EF)*e/(k*T))+1);
arvo =x*((m1)^(3/2))/pi-y*((m2)^(3/2))/pi+r-o;

*
*

% fzeron kutsu
EF = fzero(@nollattava, x0, options, E, T, EG);

H
Seuraa 
Viestejä2622
Liittynyt16.3.2005

Nopeasti katsoen vika on siinä, että käytät muuttujia x ja y funktion sisällä, mutta et välitä niitä kutsussa tai määritä niitä globaaleiksi.

Vierailija
kiltanen
kyseessä on kaava (6). Kurssin sivuilta löytyvän aputiedoston pohjalta olen rakentanut *.m-tiedoston, joka tulostaa näytölle me=9.

Mikä on tämä me? Ja mikä on siis tarkalleen ottaen ongelma: koodi kaatuu, vai ei anna oikeaa tulosta?

x=quadl(@integ_ele, EG*e, 10000,[],[],E, T, EF, EG);
y=quadl(@integ_auk, -10000, 0,[],[], E,T,EF,EG );
%******************************************
% Lauseke, jonka nollakohta etsitään
%------------------------------------------
function arvo = nollattava(E, T, EF , EG);
% Sijoita lausekkeen arvo muuttujaan arvo.
%
% Elektroni- ja aukkotiheydet sekä ionisoituneet akseptori- ja donoritiheydet
% tehtävän annon kaavasta (6) muokattuna.
r=NA*s/(4*exp((EA-EF)*e/(k*T))+1);
o=ND*s/(2*exp(-(EG-EC-EF)*e/(k*T))+1);
arvo =x*((m1)^(3/2))/pi-y*((m2)^(3/2))/pi+r-o;

*
*

% fzeron kutsu
EF = fzero(@nollattava, x0, options, E, T, EG);


Noi boldatut argumentit näyttäis ainakin menevän väärin. Mutta paha sanoa mistä homma on kiinni tältä pohjalta, pitäisi nähdä koko koodi.

Vierailija

kirjoitin funktion, jonka nollakohtaa etsitään yhteen pötköön, mutta ei auttanut. me=elektronin massa joka on määritelty aiemmin koodissa.

arvo =quadl(@integ_ele, EG*e, 10000,[],[],E, T, EG, EF)*((m1)^(3/2))/pi-quadl(@integ_auk, -10000, 0,[],[], E,T,EG,EF )*((m2)^(3/2))/pi+NA*s/(4*exp((EA-EF)*e/(k*T))+1)-ND*s/(2*exp(-(EG-EC-EF)*e/(k*T))+1));

H
Seuraa 
Viestejä2622
Liittynyt16.3.2005

Jospa tosiaan kertoisit mikä se ongalma on? Tulostaako Matlab jonkun virheilmoituksen vai onko tulos vain väärä?

Vierailija

edellisessä viestissäni oleva integraali pitäisi laskea E:n suhteen EF:n pyseässä tuntemattomana parametrina joka ratkaistaan fzerolla. tein pari muutosta koodiini, joista function arvo = nollattava(T, EG , EF); oli yksi eli tiputin E:n pois nollattavan määritelmästä koska se integroidaan pois edellisessä viestissäni olevassa lausekkeessa.

Matlab ilmoittaa minulle että koska E on tuntematon funktio tai muuttuja ei nollattava-funktiota voida käyttää fzerossa. Miten E tulisi määritellä?

integrandien määrittelyt:
%******************************************
% Integrandi elektroneille
%------------------------------------------
function arv = integ_ele(E, T, EG, EF);
arv = sqrt(E-EG*e) ./(exp((E-EF)./(k*T))+1);

%******************************************
% Integrandi aukoille
%------------------------------------------
function arv = integ_auk(E,T, EG, EF);
arv = sqrt(-E) ./ (exp((EF-E)./(k*T))+1);

H
Seuraa 
Viestejä2622
Liittynyt16.3.2005
kiltanen
Matlab ilmoittaa minulle että koska E on tuntematon funktio tai muuttuja ei nollattava-funktiota voida käyttää fzerossa. Miten E tulisi määritellä?

Tuon perusteella E ei ole määritelty siellä missä sitä tarvitaa. Matlab ilmoitti kyllä myös funktion, joten tarkasta se. Laita breakpoint funktion alkuun ja katso workspace:stä mitä muuttujia löytyy.

kiltanen
edellisessä viestissäni oleva integraali pitäisi laskea E:n suhteen EF:n pyseässä tuntemattomana parametrina joka ratkaistaan fzerolla.

Tuo ei kyllä onnistu. Kokeile toisin päin. Ratkaise ensin EF jollakin E:n arvolla ja sitten integroi => iteraatio.

Vierailija

2.3 klo 8.25 lähettämästäni viestistä näkyy funktion arvo määrittely, jonka nollakohtaa etsitään. Se mitä halutaan ratkaista on EF tietyllä lämpötilan ja siten myös EG:n arvolla, EC tunnetaan EG:n suhteen ja Ea annetaan parametrina ohjelman kutsussa.

H
Seuraa 
Viestejä2622
Liittynyt16.3.2005

Aikamoinen hässäkkä. Periaatteessa se voisi tomia. Tsekkaa se funktio mistä Matlab valittaa. Jos ei debuggaamalla vikaa löydy, yksinkertaista lauseketta, kunnes se toimii ja palauta sitten taas pala kerrallaan. Korjaa aina kun ei toimi ja lisää sitten taas pala.

Vierailija
H
Aikamoinen hässäkkä. Periaatteessa se voisi tomia. Tsekkaa se funktio mistä Matlab valittaa. Jos ei debuggaamalla vikaa löydy, yksinkertaista lauseketta, kunnes se toimii ja palauta sitten taas pala kerrallaan. Korjaa aina kun ei toimi ja lisää sitten taas pala.

Hieman samaa viestiä joudun itsekin toistamaan.

Kun kysyy apua koodaus/skriptiongelmiin, hyvä käytäntö olisi tosiaan eristää ongelmakohta mahdollisimman lyhyeen koodipätkään ja sitten kun omat taidot loppuvat, tarjota mahdollisimman lyhyt ja suoritusvalmis virheen tuottava koodi muiden tarkasteltavaksi (jos virhe ei ole ilmiselvä, maininta myös mistä on kysymys). Myös copy&paste ohjelman virheellisestä ajoyrityksesta (kometolinjalta) ja virheilmoituksesta auttaa ongelman selvittelemistä nettikommunikoinnissa!

Vierailija

Tässä olisi tämän hetkinen lauseke:
for n=1:length(T)
EF(n) = fzero(quadl( sqrt(E-EG(n)*e) ./(exp((E-EF(n))./(k*T(n)))+1), EG(n)*e, 10000,[],[], EF(n))*((m1)^(3/2))/pi-quadl(sqrt(-E) ./ (exp((EF(n)-E)./(k*T(n)))+1), -10000, 0,[],[], EF(n) )*((m2)^(3/2))/pi+NA*s/(4*exp((EA-EF(n))*e/(k*T(n)))+1)-ND*s/(2*exp(-(EG(n)-EC-EF(n))*e/(k*T(n)))+1), x0(n), options);

end

Octave herjailee tiedostostani tälläistä
`optimset' undefined near line 13 column 11
error: evaluating assignment expression near line 13, column 9
error: called from `harjoitus_1' in file `/home/kari/harjoitus_1.m
viitaten seuraavaan komentoon:

options = optimset('Display','final', 'MaxFunEvals', 1E3, 'MaxIter', 1E3, 'TolFun', 1.0D-9, 'TolX', 1.0D-9);

Miten muunnan tuon Octaven hyväksymään muotoon?

Matlab herjaa tälläistä viitaten yllä olevasta hieman poikkeavaan fzero-lausekkeeseen, sain samat herjat ylläolevastakin muodosta. nollattava on yllä kirjoitettu suoraan fzeron sisään.

??? Error using ==> fzero
FZERO cannot continue because user supplied function_handle ==> nollattava
failed with the error below.

Undefined function or variable 'E'.

Error in ==> harjoitus_1 at 16
EF = fzero(@nollattava, x0(n), options, T(n), EG(n));

Mikä oikein mättää ja miten vika korjataan?

Vierailija

Työ on saatu siihen pisteeseen, että on yritettävä keksiä miksi koodi ei toimi oikein vaan antaa aina saman tuloksen riipumatta eräiden parametrien arvoista. Lisäksi koodi ei kykene piirtämään saamaansa tulosta. Ongelma on mitä ilmeisimmin integrointirajoissa tai NA:n ja ND:n sisältävissä lausekkeissa, mutta olen kokeillut isoa määrää erilaisia integrointirajoja sen vaikuttamatta tulokseen. Siis jos NA=ND=0 koodi tuottaa oikean tuloksen, mutta kun tehtävän b-kohdassa pidetään muut parametrit samoina ja ensin NA saa arvot 10^15, 10^20, 10^25 ja ND=0 sekä seuraavaksi päinvastoin.

Ongelmallisen tästä asiasta tekee se, että koska kysesssä on harjoitustyö koodia ei ole hyvä julkaista kokonaan.

Tässä on koodin alkupää ja melkein valmiit kutsut kahteen eri tilanteeseen.

function EF = harjoitus_1(EG0,EA,EC,a,b,NA,ND,m1,m2)

e=1.60217653*10^(-19);k=1.380658*10^(-23);
h=6.6260693*10^(-34);me=9.1093897*10^(-31);
s=((((h/(2*pi))^2)/(2*me))^(3/2))*2*pi;

% harjoitus_1(1.17,0.067,0.033,0.000473,636,0,0,0.98,0.49)
% harjoitus_1(1.519,0.035,0.0058,0.000541,204,0,0,0.067,0.45)

T=5:2.5:500; EG=EG0-a.*(T.^2)./(T+b);

options = optimset('Display','final', 'MaxFunEvals', 1E3, 'MaxIter', 1E3, 'TolFun', 1.0D-9, 'TolX', 1.0D-9);

x0=EG*e/2;
for n=1:length(T)
EF(n)= fzero(@nollattava, x0(n), options, T(n),EG(n),k,e,s,m1,
m2,EA,EC, NA,ND);
end

function arvo = nollattava(T, EG , EF,k,e,s,m1,m2,EA,EC,NA,ND);

arvo = quadl(@integ_ele,EG*e,10000,[],[],T,EG,EF,k,e)*(m1^(3/2))/pi-quadl(@integ_auk,-10000, 0,[],[],T,EG,EF,e)*(m2^(3/2))/pi+NA*s/(4*exp((EA-EF)*e ./(k*T))+1)-ND*s/(2*exp(-(EG-EC-EF)*e ./(k*T))+1);

integrointirajoista vain +/-10000:tta(analyyttisissa tarkasteluissa +/- ääretön) voi muuttaa, toiset rajat tulevat tehtävänannosta.
Sitten integrandien määrittelyt jotka löytyvät tähän ketjuun 2.3 klo 18.04 lähettämästäni viestistä.

plot(T,EF)

Tämä koodi avaa piirtoikkunan muttei piirrä siihen mitään ja tulostelee ruudulle vain välejä tai herjoja joidenkin arvojen kompleksisuudesta, vaikka integrandit ovat reaalisia määrittelyalueillaan.

Miten tämän saisi toimimaan oikein?

Uusimmat

Suosituimmat