Modelleerimine COMSOL Script-iga

Kompuuterfüüsika II Praktikumid Comsol Multiphysics Skriptimine ja graafilise kasutajaliidese disainimine
  1. Käivitada COMSOL Script-i redaktor käsuga
    	edit fallingsand
    	
  2. Kustutada kõik järeltöötluskäsud reast 146 kuni 324. Järele peavad jääma kaks viimast järeltöötluskäsku.
  3. Lisaks võib kustutada %COMSOL versioonikäsud ridadelt 161 kuni 172 postplot ja postcrossplot käskude vahel.
  4. Käesolevas probleemis jääb geomeetria konstantseks, sestap pole vaja jätta geomeetriat tsükli sisse. Liigutada konstantide definitsioonid ridadelt 16 kuni 25 just fem.geom välja alla (rida 30).
  5. Redigeeri faili selliselt, et selle algus oleks kooskõlas järgnevate käskudega. Ainult esimesed read on näidatud, seega jäta mittenäidatud read puutumata. Uued ja muudetud read on kursiivis koos kommentaaridega, muudetude rea kohal on kommenteeritud originaalrida (%%-sümbolitega).
    	% COMSOL Multiphysics Model M-file
    	% Generated by COMSOL 3.5 snapshot (COMSOL 3.5.0.468,
    	$Date: 2008/08/15 16:04:11 $)
    	flclear fem
    	% COMSOL version
    	clear vrsn
    	vrsn.name = 'COMSOL 3.5';
    	vrsn.ext = ' snapshot';
    	vrsn.major = 0;
    	vrsn.build = 468;
    	vrsn.rcs = '$Name:  $';
    	vrsn.date = '$Date: 2008/08/15 16:04:11 $';
    	fem.version = vrsn;
    	% Geometry
    	g1=rect2('6e-3','14e-3','base','corner','pos',{'0','-6e-3'},
    	'rot','0');
    	g2=circ2('1e-3','base','center','pos',{'0','0'},'rot','0');
    	g3=geomcomp({g1,g2},'ns',{'g1','g2'},'sf','g1-g2','edge',
    	'none');
    	% Analyzed geometry
    	clear s
    	s.objs={g3};
    	s.name={'CO1'};
    	s.tags={'g3'};
    	fem.draw=struct('s',s);
    	fem.geom=geomcsg(fem);
    	
    %% Uus rida: iteratsiooniparameetri initsialiseerimine:

    i=1;

    %% Uus rida: FOR-tsükli algus

    for rho=2500:500:4500

    %% Konstantide väli on asetatud siia, et kaasata tera tiheduse muutused iteratsioonitsüklisse
    	% Constants
    	
    %% Muudetud rida: muuta rho_grain muutuja väärtus rho-ks
    	fem.const = {'rho_water','1000[kg/m^3]', ...
    	'eta_water','1.51e-3[Pa*s]', ...
    	'r_grain','1[mm]', ...
    	'V_grain','4/3*pi*r_grain^3', ...
    	
    %%'rho_grain','2900[kg/m3]', ...

    'rho_grain',rho', ...
    	'm_grain','V_grain*rho_grain', ...
    	'g','9.81[m/s^2]', ...
    	'Fg','-m_grain*g'};
    	% Initialize mesh
    	
    %% Muudetud read:
    %% 1. Eemaldada 'hmaxfact' (maksimaalne elemendisuuruse skaleerimisfaktor)
    %% ja 'hgrad' (elemendi kasvukiirus). Siis genereeritakse normaalne
    %% vaikevõrk, et kiirendada arvutusi.
    %% 2. Seada atomaatvõrgu parameetri väärtuseks 3 (finer)
    fem.mesh=meshinit(fem,'hauto',3);
  6. Materjaliomaduste ja ääretingimuste muutmise vajadus puudub. Iga iteratsiooni jaoks värskendab mudel automaatselt konstandi rho_grain väärtust. Pärast arvutamist tuleb lahend salvestada järeltöötluseks või muidu kustutatakse lahend automaatselt.
  7. Arvutusaja kokkuhoidmiseks seada absoluutne täpsus selle vaikeväärtusele (eemalda 'atol'-rida ja kasutada strict-tingimust solveri ajasammude jaoks.
    	% Solve problem
    	
    %% Muudetud rida:
    %% 1. Eemaldada 'atol'-rida
    > %% 2. Lisada 'tstep' parameeter väärtusega 'strict'
    	'solcomp',{'v','u','lm2','lm1','Xdot','p','X'}, ...
    	'outcomp',{'v','u','lm2','Xdot','lm1','p','X',
    	'vt','ut','lm2t','Xdott','lm1t','pt','Xt'}, ...
    	'blocksize','auto', ...
    	'tlist',[0:0.025:1], ...
    	'estrat',1, ...
    	'tout','tlist', ...
    	
    'tsteps','strict', ...
    	'linsolver','pardiso', ...
    	'uscale','none');
    	% Save current fem structure for restart purposes
    	fem0=fem;
    	
    %% Uus rida
    %% Salvestada hetkelahend järeltöötluseks
    solution(i)=fem;

    %% Kui mudelis on palju FEM-struktuure või kui mudel on üsna suur, siis pole neid
    %% mugav muutujatena salvestada. Selle asemel võib FEM-struktuuri salvestada
    %% FSLAVE-käsuga MPH-faili:
    %% fslave(sprintf('falling_sand_solution_%d',i),fem);
    %% Salvestamaks FEM-struktuuri mõnda teise kataloogi, muuta sprinft-funktsiooni
    %% vastavalt:
    %% sprintf('path_directory/falling_sand_solution_%d,i)

    %% Suurendada loendurit järgmise arvutustsükli jaoks
    i=i+1;

    %% FOR-tsükli lõpp
    end
  8. Järetöötlus ja -analüüs koosneb erinevate tulemuste kuvamises samale graafikule. Esimene alamgraafik on 2D-kiirusväli nagu saadud algselt COMSOL Multiphysics-ga. Teine alamgraafik on tera kiiruse sõltuvus ajast, mille puhul tuleb teha tsükkel üle kõigi tiheduste jaoks arvutatud lahendite. Kolmas on jõu sõltuvus tihedusest ning viimane graafik näitab stabiliseerumisaja sõltuvust tihedusest.

    %% Salvestada iteratsioonide koguarv
    max_iter=i-1

    %% Defineerida esimene alamgraafik
    subplot(2,2,1)
    	Plot solution
    	
    %% Muudetud read:
    %% 1. Eemaldada rida 'tridlim',...
    %% 2. Muuda 'title'
    %% 3. Muuda 'axis'
    	postplot(fem, ...
    	'tridata',{'U_ns','cont','internal','unit','m/s'}, ...
    	'trimap','jet(1024)', ...
    	'flowdata',{'u','v'}, ...
    	'flowcolor',[0.0,0.0,0.0], ...
    	'flowdens','velocity', ...
    	'flowlines',22, ...
    	'solnum','end', ...
    	
    %% 'title','Time=1 Surface: Velocity field [m/s] Streamline: Velocity field', ...
    	'title','Grain density 4500 kg/m3 Surface: Velocity
    	(m/s) Streamline: Velocity', ...
    	
    %% 'axis',[-0.004780328082928404,0.010780328135082468,
    -0.006700000073760748,0.008700000401586295]);
    	'axis',[-6e-3,1.2e-2,-7e-3,9e-3]);
    	
    %% Defineerida teine alamgraafik
    subplot(2,2,2)

    %% Hoida genereeritud graafikut, et võimaldada samal kuval mitut graafikut
    hold on

    %% Alustada FOR-tsüklit
    for i=1:1:max_iter

    %% Eraldada terale mõjuv kogujõud ja tera tihedus
    [Force(i),rho(i)]=postinterp(solution(i),'Fz+Fg','rho_grain', ...
    zeros(0,1),'dom',1);

    %% Defineerida RGB-kuva skaala: kasutades SIN ja COS funktsioone geneeeritakse
    %% automaatselt RGB-andmete kogum.
    R=sin(pi/2*i/max_iter);
    G=cos(pi/2*i/max_iter);
    B=sin(pi/2*i/max_iter);

    %% Muudetud rida:
    %% 1. Lisada POSTCROSSPLOT-käsu ette muutuja nimi, et see kuva seadistused
    %% oleksid ligipääsetavad
    %% 2. Muuta FEM-struktuuri definitsioon lahendimassiivis ühele liikmele
    %% 3. Muuta 'lincolor'-seadistuse parameetrit, et iga graafik oleks erinevat värvi
    %% 4. Muuta pealkiri
    %% 5. Defineerida x- ja y-pealdis
    	% Plot in cross-section or along domain
    	
    %% postcrossplot(fem,0,[1], ...
    	h=postcrossplot(solution(i),0,[1], ...
    	'pointdata','-Xdot', ...
    	
    %% 'lincolor',[0.0,0.0,0.0], ...
    	'lincolor',[R,G,B], ...
    	
    'title','-Xdot', ...
    	'title','Terminal velocity', ...
    	
    'axislabel','Time','-Xdot', ...
    	'axislabel',{'[s]','[m/s]'}, ...
    	'refine','auto');
    	
    %% Uus rida: lisada iga tiheduse jaoks legend
    set(h,'legend','on','legendstring',
    num2string(rho(i)))

    %% Täpsuse defineerimine
    tol=5e-5;

    %% Iteratsioonide maksimaalse arvu defineerimine
    Nmax=size(solution(i).sol.tlist,2);

    %% Veahinnangu algväärtustamine
    Err=1;

    %% Iteratsiooniparameetri algväärtustamine
    k=2;

    %% Tera kiiruse arvutamine iga ajasammu jaoks
    Xdot=postinterp(solution(i),'Xdot',zeros(0,1),'dom',1, ...
    'solnum','all');
  9. Stabiilse oleku saavutamiseks tuleb hinnata languskiiruse erinevust kahel järgneval ajahetkel. Kui variatsioon on väiksem kui täpsuskriteerium (või kui on jõutud viimase ajahetkeni lahendis), siis väljastab kood selle ajahetke.

    %% WHILE-tsükli algus
    while((Err>tol)&&(k<=Nmax))
    Err=(Xdot(k)-Xdot(k-1))/Xdot(k);
    t(i)=solution(i).sol.tlist(k);
    k=k+1;
    end

    %% Järeltöötluse FOR-tsükli lõpp
    end

    %% Defineerida kolmas alamkuva
    subplot(2,2,3)

    %% Kuvada jõu sõltuvus tera tihedusest
    plot(rho, Force)

    %% Lisada sellele kuvale pealkiri ja telgede pealdised
    title('Total force vs. grain density')
    xlabel('[kg/m3]')
    ylabel('[N]')

    %% Defineerida viimane alamkuva
    subplot(2,2,4)


    plot(rho, t)

    %% Lisada selle kuvale pealkiri ja telgedele pealdised
    title('Steady-state time vs. grain density')
    xlabel('[kg/m3]')
    ylabel('[s]')
  10. Skript on nüüd valmis. Salvestada see ja minna COMSOL Script-i käsureale.
  11. Veenduda, et asutakse samas kataloogis, kuhu skript salvestati. Seda saab kontrollida käsuga pwd, kataloogi muuta käsuga cd ja kataloogi sisu kuvada käskudega ls või dir.
  12. Anda käsurealt käsk fallingsand skripti käivitamiseks. Kuna kõik FEM-struktuurid salvestuvad massiivi solution, siis on võimalik kätte saada kindlale geomeetriale vastav struktuur ja re-importida see COMSOL Multiphysics-i, kus on võimalik kasutada kõiki graafilisi järeltöötlusvahendeid.

    Näiteks kolmanda raadiuse tulemuste visualiseerimiseks sisetada skripi käsurealt
    	fem3=solution(3);
    	
    COMSOL Multiphysics-i importida see kui File-Import-FEM structure: fem3.

    Kui kuva on hägune, mida võib põhjustada FEM-struktuuri importimisel mitteoptimaalset seatavad parameetrids, näiteks elemendiviimistluse (refinement) seadmine väärtusele 1, siis tuleks seada Plot Parameters-General-Element refinement: 3.


Märkus: Selleks, et FEM-struktuuri re-importimine tõrgeteta töötaks, tuleks esimese sammuna käivitada COMSOL Multiphysics , siis sealt COMSOL Script , täita skript ja re-importida struktuur juba käivasse COMSOL Multiphysics-sse.




Ülesanne 2: Esitada skripti täitmise järel saadud kuva. Esitada struktuuri re-importimisel COMSOL Multiphysics-i saadud kuva.



This document was translated from LATEX by HEVEA.

Licensed under the GNU Free Documentation License