Incalzirea Materialelor In Camp de Microunde

Cuprins

CAPITOLUL I INCALZIREA MATERIALELOR IN CAMP DE MICROUNDE

I.1 Introducere

I.2 Transferul de caldura

I.3 Specificul încălzirii în câmp de microunde

I.4 Fenomenul de difuzie termică

CAPITOLUL II PROGRAMUL DE MODELARE CU ELEMENT FINIT GFEM

II.1 Introducere

II.2 Mod de organizare

II.3 Utilizare program

II.4 Studii de caz

CAPITOLUL III MODELAREA FENOMENELOR

DE CAMP ELECTROMAGNETIC

III.1 Introducere

III.2 Metode integrale și descompunere modală

III.3 Metoda elementului finit

III.4 Modelul de cavitate cu corp de proba

III.5 Determinarea profilelor de camp electromagnetic

CAPITOLUL IV MODELAREA FENOMENELOR DE CAMP TERMIC

IV.1 Introducere

IV.2 Determinarea profilelor de camp termic

IV.3 Programul informatic MacGFEM

IV.4 Concluzii

CAPITOLUL V SOFTWARE UTILIZAT ÎN CERCETAREA DIN DOMENIUL MICROUNDELOR

V.1 Instrumente pentru modelarea câmpului magnetic, electric și termic

V.2 Analiza statistica a datelor experimentale

V.3 Software pentru studiul meta-analitic al literaturii de cercetare

V.4 Resurse disponibile în domeniu

BIBLIOGRAFIE

ANEXE

Cuvânt înainte

O cercetare solidă și un studiu aprofundat al fenomenelor electromagnetice, are nevoie de modele matematice , care pe baza unor rezultate și masurători practice , pot recrea virtual realitatea .

Modelarea și simularea fenomenelor de câmp electromagnetic și termic ale corpurilor aflate în câmp de microunde, oferă posibilitatea unor investigații știintifice rapide și suficient de precise, astfel încât cercetătorul are posibilitatea permanentă de a valida rezultatele experimentale cu cele teoretice.

Ingineria și informatica au creat acest instrument minunat , care este modelarea și care reprezintă un pas important în cercetarea integrată.

Autorii

CAPITOLUL I Incalzirea Materialelor in Camp de Microunde

Cuprins

I.1 Introducere

I.2 Transferul de caldura

I.3 Specificul încălzirii în câmp de microunde

I.4 Fenomenul de difuzie termică

I.1 Introducere

Fenomenele fizice care însoțesc transferul de energie electromagnetică într-un material, sunt complexe iar procedeele de abordare rămân dificile, deoarece sunt legate de numeroși factori interdependenți și variabili în cursul derulării procesului.

În acest capitol se va aborda modul de incalzire al corpurilor dielectrice aflate in camp de microunde, transferul de căldură și ecuațiile care guvernează acest proces.

Legile electromagnetismului sunt formulate cu ajutorul ecuațiilor lui Maxwell, care exprimă relațiile între sarcini electrice, curenți și câmpuri electromagnetice:

(I-1)

unde: – câmpul electric, – inducția electrică, – inducția magnetică, – câmpul magnetic.

– reprezintă densitatea volumică de curent

– reprezintă densitatea de sarcină

Răspunsul mediului de propagare , la excitația prin câmpul electromagnetic, se traduce prin relația care cuplează vectorii , ; , și ( , ) și care rezultă din proprietățile electrice și magnetice ale mediului.

(I-2)

unde – conductivitate electrică

– permitivtate electrică

– permeabilitate magnetică

= r0 – unde r -permitivitate relativă a mediului

0 -permitivitatea vidului [ F/m]

= r0 -În mediile magnetice =0

– unde 0 permeabilitatea magnetică a vidului

Pentru încălzirea cu microunde ca medii magnetice se pot folosi feritele.

Rescriind ecuațiile lui Maxwell, în afara sursei de radiații și ținând cont de mediu, obținem în complex:

(I-3)

Proprietățile disipative ale unui material sunt determinate de coeficientul și de partea complexă a lui r . Acești coeficienți caracterizează puterea absorbantă a unui material aflat în câmp de microunde și prin aceasta capacitatea lui de a se încălzi.

Sarcinile electrice libere, în prezența unui câmp electromagnetic dau naștere la un curent de conducție de forma :

(I-4)

Oscilațiile acestor sarcini sunt împiedicate de existența particulelor fixe și datorită ciocnirilor care au loc se produce încălzirea conform Legii lui Ohm.

Pe de alta parte, în prezența unui câmp exterior , moleculele tind sa se orienteze în așa fel încât momentul lor dipolar să fie paralel cu

. Datorita însă frecărilor cu celelalte particule, apare o întârziere a mișcării moleculelor, rezultând un defazaj între și care corespunde unei permitivități complexe.

= ’-j” (I-5)

Acest mecanism se numește relaxare dielectrică, care se traduce prin prezența unui curent fictiv:

(I-6)

Luînd în considerare aceste elemente și revenind la ecuațiile lui Maxwell obținem :

Termenul ” este echivalent cu o conductivitate.

Se definește unghiul de pierderi prin: iar cînd este neglijabil, avem:

În tabelul de mai jos se prezintă câteva valori ale lui r și tg pentru diferite materiale, [11]

Tab.I.1. Constante de material

Se poate concluziona că sub acțiunea câmpului electric, mișcările sarcinilor libere și ale moleculelor polare, provoacă ciocniri și frecări cu particulele care le înconjoară, creând o încălzire a mediului. Acesta este principiul încălzirii în câmp de microunde. [19], [20].

Încălzirea în câmp de microunde, poate avea de asemenea la origine pierderile magnetice ale mediului de propagare. În anumite medii, cum ar fi feritele, permeabilitatea este complexă =’-j” , câmpul magnetic contribuind la disiparea energiei microundelor.

I.2 Transferul de caldura

I.1.1. Energia termică și transferul de căldură

Se poate defini ca încălzire clasică, modelul de încălzire indirectă, în care energia este transmisă de la sursa de căldură, cum ar fi : rezistori, becuri infraroșii, către corpul ce trebuie încălzit, respectând legile termodinamicii clasice.

În acest caz, corpul primește energia, numai prin suprafața sa, lucru care nu se întâmplă în cazul încălzirii în câmp de microunde.

Energia termică, sau cantitatea de căldură, este suma energiilor mecanice ale particulelor, datorată agitației moleculare.

Energia primită de un volum v al corpului, are ca efect creșterea agitației particulelor și deci creșterea temperaturii, sau poate contribui la schimbarea stării de agregare, ori la o reacție chimică [138 ]. Încălzirea corpului cu T, este proporțională cu cantitatea de căldură primită Q :

Q = m c T ( I – 7 )

unde c – căldura specifică, m – masa corpului de probă.

Când are loc o schimbare de stare sau o reacție chimică endotermică, o parte a acestei energii este utilizată pentru a produce reacția:

Q = m L ( I – 8)

unde L – căldură latentă

Transmiterea energiei către corp, are loc fie prin transfer de energie termică, prin conducție și convecție, fie prin intermediul radiației în infraroșu.

Conducția

La contactul dintre două solide diferite, sau în interiorul aceluiași solid, agitația moleculelor se transmite din aproape în aproape, prin ciocniri cu moleculele vecine.

Moleculele din straturile superficiale primesc energie și o transmit prin conducție către straturile din adâncime, acesta fiind un fenomen lent, care tinde spre echilibrul static, definit de temperatura T.

Viteza de difuzie a căldurii, depinde de conductivitatea termică, definită prin constanta k. Aceasta are valori diferitem de la 395 W/m – pentru cupru, la 0,04 W/m– pentru vata de sticlă, care este un bun izolant termic.

Reacția în lanț de la o moleculă la alta, este mult accentuată de mișcarea electronilor liberi, ceea ce explică proprietățile bune conductoare, pe care le au metalele.

În cazul dielectricilor, pe motivul slabei conductivități termice a acestora, există dificultăți în încălzirea lor cu metode clasice, motiv pentru care încălzirea lor în câmp de microunde prezintă interes.

Analiza conducției termice este bazată pe legea lui Fourier, de unde rezultă ecuația căldurii, iar transferul de căldură se face de la punctul mai cald spre punctul mai rece, fluxul de căldură fiind caracterizat de :

Convecția

Convecția reprezintă transferul de căldură, care are loc în prezența unui fluid, prin ciocnirile dintre particulele fluidului și dintre fluid și corp, sau ca o combinație între cele două.

Studiul convecției este legat de curgerea fluidelor și vom reține doar expresia fluxului termic care traversează peretele unui solid :

Acesta depinde de mulți factori, cum ar fi, viteza fluidului, caracteristici ale fluidului, tipul de curgere, pereții [10].

Cazul descris este întâlnit într-un cuptor cu microunde, în care o sarcină este înconjurată de aer și Hc poate varia între 5- l50 W/m2 C , în funcție după cum aerul este liniștit sau în mișcare.

Radiația în infraroșu

Un corp încălzit, emite prin suprafața sa o radiație electromagnetică în domeniul infraroșu ( 0,76 – 100 m), deasemenea un corp aflat într-un mediu traversat de radiații infraroșii se încălzește.

Explicația ține de faptul că trecerea atomului dintr-o stare energetică în alta, se face prin absorbție sau emisie de radiații. Relația dintre energia primită sau produsă și frecvența de radiație a fost stabilită de Plank [13],

Trebuie remarcat că energia fotonului, corespunzând radiației în infraroșu, este mică în comparație cu radiația ultravioletă, radiația x sau gama și nu poate declanșa procese de alterare a nucleului atomilor, sau efecte fotochimice. În cazul microundelor energia fotonului este mai mică.

În ceea ce privește adâncimea de pătrundere a unei radiații infraroșii, aceasta este de numai câțiva milimetri. Prin urmare numai suprafața mediului va avea rolul de receptor sau emițător de infraroșii.

Legile radiației sunt studiate, pornind de la un mediu ipotetic, corp negru, care este un emițător și un absorbant ideal, adică factorul său de absorbție este egal cu 1, indiferent de temperatura și frecvența radiației incidente , distribuția spectrală a emisiei depinde numai de temperatura corpului .

Legile emisiei, bazate pe statistica Bose-Einstein, exprimă relațiile dintre trei parametrii: putere radiată, temperatura corpului și lungimea de undă a radiației.

În figura I.1, se poate observa, la diferite temperaturi cum variază puterea radiată, precum și existența unui maxim al puterii pentru o valoare a lungimii de undă M. [30 ]

Fig. I.1. Puterea radiată în funcție de lungimea de undă

S-a constatat că radiațiile în banda de frecvență a microundelor sunt emise la temperaturi apropiate de cele ale corpului uman. De aici, s-a ajuns la o aplicație medicală a microundelor, anume de a detecta anumite leziuni interne ale țesuturilor, urmărind variațiile acestor radiații.

Legea Stefan-Bolzmann fiind de forma :

unde E – factor de emisivitate al corpului , inferior lui 1, depinzând de starea suprafeței corpului și de temperatura sa.

Pentru fiecare direcție de propagare și lungime de undă, factorul de emisie este egal cu factorul de absorbție, (conform Legii lui Kirchoff), dar puterea de absorbție totală nu este echivalentă cu puterea de emisie totală. Explicația este că factorul de absorbție al unui corp depinde și de natura radiației primite (natura și caracteristicile obiectului emitor).

Considerând acum două medii cu temperaturile T1 pe suprafața emitoare și T2 pe suprafața receptoare, generalizând legea lui Stefan-Bolzmann, fluxul de căldură transmis, este:

Pentru temperaturi ridicate, peste 800 C, densitățile de putere schimbate prin radiații infraroșii sunt mult mai importante decât cele schimbate prin convecție. Dimpotrivă, la temperaturui mai mici, 200 a temperaturi apropiate de cele ale corpului uman. De aici, s-a ajuns la o aplicație medicală a microundelor, anume de a detecta anumite leziuni interne ale țesuturilor, urmărind variațiile acestor radiații.

Legea Stefan-Bolzmann fiind de forma :

unde E – factor de emisivitate al corpului , inferior lui 1, depinzând de starea suprafeței corpului și de temperatura sa.

Pentru fiecare direcție de propagare și lungime de undă, factorul de emisie este egal cu factorul de absorbție, (conform Legii lui Kirchoff), dar puterea de absorbție totală nu este echivalentă cu puterea de emisie totală. Explicația este că factorul de absorbție al unui corp depinde și de natura radiației primite (natura și caracteristicile obiectului emitor).

Considerând acum două medii cu temperaturile T1 pe suprafața emitoare și T2 pe suprafața receptoare, generalizând legea lui Stefan-Bolzmann, fluxul de căldură transmis, este:

Pentru temperaturi ridicate, peste 800 C, densitățile de putere schimbate prin radiații infraroșii sunt mult mai importante decât cele schimbate prin convecție. Dimpotrivă, la temperaturui mai mici, 200 C, aceste densități sunt comparabile [34].

Considerând condițiile asupra fluxului de căldură, identice, pentru diferite moduri de transfer, prin convecție și prin radiație, se poate scrie :

unde H- coeficientul global de transfer, determinat experimental.

Se poate concluziona că, oricare ar fi modul de transfer de căldură, prin conducție, convecție sau radiație, sursele de căldură fiind întotdeauna situate pe periferia produsului, avem de a face deci cu o încălzire de suprafață.

În încălzirea în câmp de microunde, cu toate că energia este transmisă direct în interiorul produsului, nu pot fi neglijate condițiile termice pe suprafață, unde pot totuși exista surse de căldură. În această situație trebuie să ținem cont și de faptul că :

încălzirea este mixtă – există o încălzire complementară (infraroșu, aer cald), care se adaugă încălzirii cu microunde.

există pierderi de energie la suprafață- fie printr-o izolație termică insuficientă, fie datorită radiației produsă de corpul încălzit, fie datorită evaporării (cazul apei), care evaporându-se pe suprafață, disipă energie.

Se poate aprecia deci că încălzirea clasică , care are loc de la periferia corpului, spre interior, este o încălzire lentă dar omogenă.

I.3 Specificul încălzirii în câmp de microunde

Făcând o comparație între încălzirea unui corp cilindric, aflat într-un mediu care îl încălzește și apoi aflat într-un câmp de microunde, ne putem da seama că există un complex de fenomene care rezultă din transferul de energie prin microunde, așa cum se vede în Fig. I.2.

Figura I.2. Dispersia căldurii

-Incălzirea este rapidă.

Explicația rezidă în faptul că straturile din adâncime ale corpului captează direct căldura și cum energia microundelor este disponibilă imediat, nu există timpii morți din încălzirea clasică, necesari pentru ajungerea la o anumită temperatură a cuptorului.

Dielectricii sunt în majoritatea cazurilor, slabi conducători de căldură, durata încălzirii cu metoda clasică fiind mare, în câmp de microunde însă, captează bine energia electromagnetică și datorită specificului acestei încălziri, se pretează foarte bine pentru aceste medii.

– Incălzirea este neuniformă.

Datorită faptului că mediile sunt eterogene și limitate, încălzirea cu microunde nu poate fi uniformă, fenomenele de propagare ducând la mari diferențe în ceea ce privește repartiția puterii absorbite.

Există metode de îmbunătățire a termodinamicii în ansamblu: fie se încearcă o variație de la un moment la altul a distribuției de putere , fie excitația de putere se face intermitent, lăsând timp pentru difuzia termică, fie prin asocierea unei alte surse de energie ( aer cald).

– Starea suprafeței corpului.

Temperatura în interiorul corpului este mai mare decît pe suprafața lui, corpul încălzind cuptorul (prin radiațe sau prin convecție). Pentru a putea încălzi suprafața corpului, aceasta trebuie izolată termic. În stratul exterior al corpului are loc un fenomen de disipare termică suplimentară, datorită evaporării apei conținute. În lipsa unui curent de aer uscat, evaporarea apei se face cu dificultate.

– Ambalarea termică.

Ambalarea termică se produce atunci când apare un punct cald și când factorul de pierdere al mediului crește cu temperatura. Rapiditatea încălzirii duce la un dezechilibru termic, accentuat de fenomenul de ambalare termică [15]. Acest proces este foarte sensibil în cazul decongelării sau al unor reacții chimice.

– Selectivitatea.

Încălzirea diferențiată obținută în câmp de microunde, se poate folosi la obținerea unor tratamente selective.

Ca exemple pot fi date, lipirea lemnului unde numai cleiul care prezintă pierderi dielectrice, va fi încălzit și polimerizat ; anumite sterilizări în care paraziții ce trebuie eliminați, captează energia microundelor etc. Deasemenea un deosebit interes prezintă selectivitatea în reacțiile chimice în medii eterogene, în situația în care catalizatorul este disipativ, energia este absorbită de acesta și reacția va fi cu atât mai accelerată.

I.4 Fenomenul de difuzie termică

Un material de volum v, este aflat într-un câmp electromagnetic și dorim să-i aflăm profilul de temperatură T1 (M,t). Volumul este plasat într-un alt mediu ale cărui caracteristici de temperatură, presiune și umiditate sunt variabile în timp. Temperatura T1 (M,t), este soluția

ecuațiilor căldurii care exprimă transferul

prin conducție a energiei termice în interiorul

volumului v și care trebuie să satisfacă

condițiile la limită pe frontiere.

Studiul conducției termice se bazează pe Legea lui Fourier, care exprimă în fiecare punct al mediului, relația dintre fluxul de căldură, prin unitatea de suprafață, J și variația temperaturii.

unde K- conductivitatea termică a mediului.

Semnul “- “ semnifică faptul că transferul de căldură se face în sensul temperaturilor descrescătoare, de la regiunile calde spre cele reci, în scopul egalizării temperaturilor.

I.3.1. Ecuația căldurii

Într-un punct oarecare al materialului încălzit, conservarea energiei este descrisă de următoarea formulă :

Rezolvarea acestei ecuații permite calculul distribuției temperaturilor în produsul tratat la un moment dat, precum și viteza de variație a temperaturii [73].

Dacă k este invariant, atunci :

I.3.2. Termenul sursă

Termenul P(M,t), reprezintă puterea furnizată sau pierdută simultan de către sursele interne de căldură existente în mediu. Dacă P este pozitiv, atunci exprimă puterea furnizată mediului, dacă P este negativ, atunci reprezintă puterea pierdută.

Puterea furnizată materialului este de fapt puterea disipată la momentul t în puncutl M. Pe unitate de volum vom avea:

Acest termen este variabil în volumul produsului încălzit și în timpul tratamentului. Într-adevăr distribuțuia câmpului este neuniformă, datorită prezenței undelor staționare și eterogeneității mediului, care se transformă în timpul încălzirii. Termenul sursă P, este deasemenea dependent de aplicator.

Vom putea vedea că într-o incintă închisă, câmpurile sunt impuse de geometria aplicatorului dar mai ales că există un cuplaj între sarcină și aplicator.

În continuare abordăm două exemple privind forma de variație a câmpului

electric, într-unb ghid de undă și într-un cuptor casnic, obținute prin utilizarea metodei elementului finit [56],

Fig II.5. Ghid de undă cu sarcină

Fig. II.6. Aplicator dreptunghiular cu sarcină

Fig. II.7. Ghid vid. Distribuția lui E

Fig. II.8. Ghid cu eșantion. Distribuția lui E

Figura II.9.a. Influența eșantionului asupra lui Ez

Figura II.9.b Influența eșantionului asupra lui Ex

Fig. II.10. Repartiția puterii în dielectric

I.3.3. Transferul termic la periferia materialului

Studiul condițiilor la limită asupra temperaturii, revine la cercetarea condițiilor asupra surselor de căldură externe la suprafața materialului.

În funcție de natura produsului încălzit și de mediu, condițiile la limită pe suprafața s , delimitând un volum v, sunt diferite .

A. Se impune o valoare a temperaturii la suprafață

TS=T0 (I- 23)

În acest caz, schimburile termice între cele două medii sunt efectuate în ambele sensuri , dar sunt echivalente.

Schimb de energie cu mediul înconjurător

Schimburile de căldură se efectuează într-un singur sens de la mediul mai cald spre mediul mai rece. Cantitatea de căldură schimbată cu exteriorul prin suprafață, este :

unde H- caracterizează materialul și starea suprafeței sale.

Determinarea teoretică a temperaturilor într-un material încălzit cu microunde se face prin rezolvarea ecuației încălzirii, ținând cont de condițiile la limită pe pereți. Termenul sursă va trebui calculat prin rezolvarea ecuațiilor de undă.

CAPITOLUL II Programul de Modelare cu Element Finit GFEM

Cuprins

II.1 Introducere

II.2 Mod de organizare

II.3 Utilizare program

II.4 Studii de caz

II.1 Introducere

Analiza numerică pentru ecuațiile diferențiale parțiale reprezintă 90% din calculul numeric științific de azi. Majoritatea cazurilor pot fi soluționate de către un program cu Element Finit. Există două tipuri de instrumente de element finit:

– instrumente simple, care sunt de ajutor, dar necesită multă pregătire de date, precum triangulații.

– instrumente puternice profesionale care sunt deseori costisitoare și pentru care timpul învățării este mult prea mare.

Astfel, Gfem dorește să fie la mijloc: ieftin, puternic, cu un scurt timp de învățare.

II.2 Mod de organizare

Gfem este scris în jurul unui nucleu numit FreeFEM. Acesta citește un program, îl execută și produce fișiere și grafice / diagrame , nu este interactiv; dar Gfem este! Este coordonat de meniu și de mouse. De exemplu, se poate scrie un program din cadrul Gfem și să se execute doar parte din el, sau să fie modificat și executat fără a părăsi aplicația. Graficele / diagramele pot fi mărite, salvate, tipărite. Cum se aude, FreeFEM este gratuit ( cu bunăvoința Numerica ) și portabil pe PC – uri și pe Unix-X11-uri astfel să vă puteți executa programul dezvoltat cu Gfem și pe aceste aparate; de asemenea, poate fi dat de asemenea și studenților .dar Gfem este un produs comercial complet interactiv ( „g” vine de la Grafice și Limbaj ) .

Fig.1 Cele 3 ferestre ale Gfem pe un Macintosh . Fluxul

potențial în jurul unei configurații duble NACA.

Fundamental, Gfem lucrează pe funcții scalare de 2 variabile x, y, definite pe un domeniu . Domeniile sunt reprezentate prin triangulații, iar funcțiile prin aproximații liniare continue pe porțiuni pe triangulații. Operațiile majore la funcții sunt formulele algebrice, derivate, convecții, inversare a operatorilor eliptici.

De exemplu diagrama(dx(sin(x+y))) va afișa nivelul curbelor interpolării lineare ale derivatei x ale sin(x+y) pe un domeniu definit mai devreme de către o triangulație.

Triangulația

Generarea automată a mesh-ului folosind elemente triangulare este efectuată pentru un domeniu a cărui limită / frontieră este introdusă analitic.

De exemplu, următoarea comandă va triangula circumferința unitară cu o gaură în interiorul radius-ului 0.1 și un maxim de 820 vertice. Numărul de referință al frontierei exterioare este egal cu 1; este discretizat cu 80 puncte. Frontiera interioară are un număr de referință 2 și 20 puncte. Finețea mesh-ului este controlată de mărimea celei mai apropiate limite marginale.

deuxpi:=8*atan(1)

limita (1,0,deuxpi,80)

(x:=cos(t); y:=sin(t) );

limita (2,0,deuxpi,20)

(x:=0.3+0.5*cos(-t); y:=0.5*sin(-t) );

Buildmesh (820);

Fig.2 Triangulare, exemplu

La fel pot fi tratate și domeniile cu colțuri. Domenii foarte complexe pot fi generate cu exactitate, așa cum este arătat pe imaginile de mai jos.

Operațiunile simple

Funcțiile sunt generate prin atribuire precum în f=x*y; funcțiile pot fi operate cu ’*,/,+,-,’, funcțiile obișnuite sin, cos … și operatori noi precum ’plot’, ’dx’… De exemplu:

F=sin(x)*y; g=dx(f); plot(tan(x)*g-f);

Condiții de frontieră pentru EDP-uri ( ecuații cu derivate parțiale)

Luați în calcul cazul unui flux potențial asupra unui cilindru, rezolvat cu o funcție de curent v. referindu-ne la triangulația de mai sus, vom scrie:

Onbdy(1) v = u2*y-ul * x;

Onbdy(2) v = 0;

Neumann și robin / la fel sunt tratate și condițiile de frontieră Fourier

Rezolvatorul

În forma sa prezentă Gfem are un rezolvator general pentru ecuațiile cu derivate parțiale eliptice scalare:

Aceasta folosește metoda elementului finit de gradul unu pentru a construi un sistem linear și apoi o factorizare Gauss pentru a o rezolva. Acestea sunt numite prin scris

solve(v)

id(v)*b+dx(v)*a1+dy(v)*a2-laplace(v)*nu=f

în timp ce matricele sunt salvate, este posibil să rezolvați mai multe ecuațiile cu derivate parțiale dintr-o dată, fără refactorizare.

Puterea lui Gfem este pe probleme dependente de timp pentru că este posibil să facă un ciclu. Aici sunt pași de 10 timpi pentru ecuația căldurii:

dt:=0.1; v = 1;

onbdy (1) v = 0;

solve ( v, 1 )

id (v ) / dt – laplace ( v ) * nu = v / dt;

iter (10)

(onbdy (1) v=0

solve ( v, -1 )

id (v ) / dt – laplace ( v ) * nu = v / dt;)

2.Grafice

Sunt disponibile următoarele configurații:

Afișaj de trei ferestre ajustabile ca mărime care conțin:

– programul

– triangulația

– nivelul curbelor sau harta color sau secțiune transversală 1D

Ferestrele grafice pot fi mărite fără limită. Conținutul fiecărei ferestre poate fi stocat în format text. Toate graficele pot fi copiate și mutate într-o altă aplicație ca și fișiere EPFS pict ( Macintosh ). Toate graficele pot fi, de asemenea depozitate în fișiere Postcript cu antet pentru a le include în fișiere normale TEX într-o formă rezonabilă.

3. Anexă: Caracteristici ale FreeFEM

FreeFEM este un program pentru rezolvarea Ecuațiilor cu derivate parțiale ( EDP). El funcționează pe Macintosh, PC și Unix-X11.

Există trei pași pentru soluționarea unei EDP

Pregătirea datelor unde frontierele domeniului al EDP sunt definite împreună cu condițiile de frontieră și alte date.

Pasul de execuție care implică construirea triangulației ,construcția sistemelor lineare sau non-lineare și soluțiile lor

Afișarea datelor

FreeFEM are un limbaj pentru pasul 1 care folosește cuvinte cheie pentru a lansa pasul 2 și 3.

solve (u) rezolvă EDP și se întoarce u

plot (u) afișează liniile de nivel ale u

save (u, ’file.txt’ ) scrie u pe fișierul de disc ’file.txt’ în format text

Celelalte cuvinte cheie ale limbajului sunt:

„sin”, „cos”, „atan”, „exp”, „log”, „abs”, „sqrt”, „min”, „max”.

„if”, „then”, „else”

„iter”, „border”, „buildmesh”, „onbdy”, „dnu”, „laplacien”, „id”, „dx”, „dy”, „load”, „loadmesg”, „plot”, „plot3d”, „execute”, „saveall”, „savemesh”.

În fine, Gfem are 6 variable rezervate:

„x”, „y”, „t”, „ib”, „region”

și folosește următoarele delimitatoare: „ ( () ) , ; . ” și următorii operatori:

,”*”, „/”, „+”, „-„, „ „, „”, „”, „”, „”, „==”, „and”, „or”, „=”, „:=”.

Sintaxa limbajului este, în linii mari, ca și cea a limbajului Pascal.

FreeFEM citește programul său introdus dintr-un fișier de pe disc. Pe Unix sau PC – DOS numele programului este primul argument, ce apare după numele „freefem”. Pe Macintosh apare un fișier de dialog odată ce programul este lansat.

Lansați programul. Între fiecare imagine, APĂSAȚI MOUS-ul.

4. Modul de lucru

– Frontierele sunt definite prin părți, fiecare fiind dată de ecuația sa x(t), y(t), și etichetată cu un număr ib. Luând un număr finit de valori N, parametrul t definește N puncte de frontieră care vor fi introduse în triangulator.

– triangulatorul folosește metoda Dealunay – Voronoi. Până ce nu este atins numărul maxim specificat de puncte, el generează vertice interioare / interne pe mijlocul marginilor până ce toate marginile sunt mai scurte decât media greutăților punctelor finale ale sale. Greutatea fiecărui vertex de frontieră este lungimea medie ale celor două margini de frontieră din jurul său. Greutatea atribuită unui nou vertex interior este media greutății celor două puncte capăt ale marginii a cărei mijloc este. Punctele mesh sunt re-numărate de către un algoritm revers Cuthill – McKee.

– Rezolvatorul EDP este o Metodă de Element Finit (MEF ) de ordinul 1 pentru o ecuație scalară cu două tipuri de condiții de frontieră, Dirichlet sau Robin.

Sistemul linear este rezolvat prin factorizarea benzii Gauss. Toate formulele cvadraturice presupun ca datele sunt pe porțiuni și continue.

Mod de lucru cu interpretatorul

– precum majoritatea interpretatorilor, are o parte lexicală și una sintetică. În partea lexicală, sursa este ruptă în atomi lexicali și recunoscută ca simboluri ( vezi simbolul enum în fișierul sursă lexical.h ). De exemplu, dacă primul caracter este o cifră atunci acesta este un număr iar tipul de simbol asociat este „cste”. Acest lucru este efectuat de către funcția „nextsym”. În plus, construiește un tabel de constante și variabile ( care pentru utilitate conține și cuvintele rezervate ale limbajului ).

– Analizatorul lexical este o funcție numită de către analizatorul de sintaxă. Prin urmare, a doua este rutina principală în program, cu excepția câtorva inițializări; numele său este „instrucție”.

Analiza sintaxei este acționată de către regulile de sintaxă datorită faptului că limbajul este LL(1). Astfel există funcția C pentru fiecare non – terminal. – programul nu generează un cod de obiect dar în schimb, generează un graf orientat aciclic . De exemplu, derivarea x*y ar genera un graf orientat aciclic cu rădăcina * și două ramuri „x” și „y”. Grafii orientați aciclic sunt C-struct cu patru indicatori de adresă , pentru cele 4 ramuri ( aici două ar indica NULL ) și un simbol pentru operator. Funcția C care construiește grafi orientați aciclic, este numită plante. La final, programul este transformat într-un graf complet, și pentru a executa programul, trebuie să faceți doar un singur lucru: evaluați operatorul la rădăcina grafului.

– Programul de evaluare este efectuat de către funcția C „eval”. Acesta se uită la simbolul rădăcinii și execută operațiunea corespunzătoare. Aici ar face

întoarcere „valoare indicată de L1”*”valoare indicată de L2”

dacă L1 și L2 ar fi adresele celor două ramuri.

– Tipul jocului este de a asocia un graf fiecărei operațiuni. De exemplu, când valoarea variabilei „x” este necesară, acest lucru este de asemenea făcut de către un graf care are operatorul „oldvar” ca rădăcină. Cel mai mare artificiu dintre toate este instrucțiunea compusă „…;…;”. Aici „” este considerat ca un operator cu o ramură pe instrucțiunea curentă și cu o ramură pe următoarea. Similar pentru instrucțiunea if … then … else.

Să presupunem că cineva dorește să adauge o instrucțiune în FreeFEM, iata ceea ce trebuie făcut:

– Asigurați-vă că sintaxa este LL(1) și nu intră în conflict cu cea veche. Adăugați cuvinte rezervate structurii de date cu „installe”. Datorită faptului că se va adăuga un nou simbol, actualizați lista de simboluri ( o structura enum).

Adăugați funcțiile C pentru analiza sintaxei în conformitate cu diagramele. Modificați eval-ul adăugând la „switch” noul caz .

II.3 Utilizare program

Scopul

Gfem îndeplinește trei sarcini:

Creează, editează și afișează o triangulație

Creează și afișează funcții pe o triangulație

Rezolvă Ecuații diferențiale Parțiale

Gfem este rulat de programe citite de pe disc sau scrise în ferestre de editare ( edit windows ).

Într-un mediu integrat precum MacGfem, PCGfem sau XGfem, pentru executarea unui program, puteți folosi meniul Run sau Select parte din text și apăsați ENTER ( nu return ). Fără mediul meniu – fereastră, numele fișierului programului este citit în momentul lansării.

Acest program simplu: ( Frontieră = border )

border ( 1,0,6.28,20 ) ( x:=cos(t); y:=sin(t) );

buildmesh ( 300 );

f=x*y; plot (f);

va crea o triangulație pentru cercul unitate cu 20 puncte pe frontiera care a primit eticheta 1, și un maxim de 300 vertice, apoi afișează liniile de nivel ale funcției f=x*y;

Dacă programului anterior este adăugat:

onbody (1)u = 0; solve (u) – laplace (u) = f

Atunci este generată soluția

– ∆ u = f în ,u∂ = 0

Convenții , elemente de condiționare

O triangulație, T, este un set de triunghiuri, fiecare având 3 vertice, acoperind un domeniu referit ulterior ca . Fiecare triunghi și fiecare vertex are un număr și în plus triunghiurile au o regiune de număr de referință a regiunii în timp ce verticele au de asemenea un număr de referință frontieră, ib, care Identifică Frontierele.

Fiecare vertex al triangulației are un număr numit iv.

O funcție matrice este dispunere de valori pe verticele T-ului. Reprezintă o funcție continuă lineară pe porțiuni pe .

Într-un mediu integrat Gfem are noțiunea Proiectului: fișiere arătate pe desktop ca și documente cu o icoană precum un dreptunghi triangulat.

Un proiect conține un program, ultima triangulație creată și o funcție array, adică ultima funcție afișată.

Triangulații

Pentru a crea o triangulație trebuie sa procedam astfel :

Să deschidem un proiect

Să citim o triangulație veche salvată în format text

Să executam un program care să conțină keyword = cuvânt cheie buildmesh.

Să creați, desenînd frontiera – ului și activați meniul „Triangulate” ( doar Macintosh ).

În medii integrate, odată create, triangulațiile pot fi afișate, depozitate pe disc în format Gfem sau în format text sau chiar o mărire ( zoom ) a afișajului sau grafic poate fi salvat în format Postscript și inclus într-un fișier TeX.

Gfem depozitează pentru fiecare vertex numărul său și numărul de referință al frontierelor sale ( =0 dacă nu e pe o frontieră ) și pentru fiecare triunghi numărul său și numărul de referință regiune. Numărul marginilor nu este depozitat, dar pot fi recuperate de program dacă e necesar.

Funcții și scalare

Funcțiile sunt fie citite sau create.

Funcțiile pot fi citite dintr-un fișier, dacă valorile sale la verticele triangulației sunt depozitate în format text. ( deschideți un exemplu cu un editor text pentru a vedea formatul ).

Funcțiile pot fi create prin executarea unui program. Instrucții precum f=x*y înseamnă de fapt că f(x,y) = x*y pentru toate coordonatele x și y ale verticelor. Aici x și y se referă la coordonatele din domeniul reprezentat de către triangulație.

Cele mai utilizate funcții ce pot fi folosite:

max, min, abs, plus, minus,

cos, sin, tan, acos, asin, atan,

cosh, sinh, tanh, log, exp, sqrt

funcțiile pot fi create cu alte funcții anterior definite precum în g=sin(x*y); f=exp(g).

Operatori:

and, or, , , , , , +, -, *, /, ^

Două alte variabile pot fi folosite pe lângă x și y: ib și region ( vezi variabilele geometrice ).

Funcțiile create de un program sunt afișate doar dacă este folosit cuvântul cheie plot ( ).

Derivatele funcțiilor pot fi create prin cuvinte cheie dx și dy. Dar ele sunt interpolate astfel că rezultatul este linear continuu pe porțiuni.

Similar, operatorul de convecție convect ( f, u, v, dt ) definește o nouă funcție care este aproximativ f (x – u (x,y)dt, y – v ( x,y)dt). Luați aminte că declarația f = convect (f, u, v, dt ) ar da un rezultat incorect datorită faptului că modifică f la unele vertexuri unde vechea valoare este încă necesară pe mai târziu. Este necesar să scrieți:

G = convect ( f, u, v, dt ); f = g;

Scalarele sunt de asemenea importante pentru a crea funcții. Din moment ce nici o funcție array nu este atașată la un scalar, simbolul := este de folos pentru a le crea, ca și în pi:= 4*atan(1); f = x*cos (pi*y); aici f este o funcție dar pi este un scalar.

Este posibil a se evalua o funcție la un punct ca și în a:=f(1,0). Aici valoarea lui a va fi 1 pentru că f(1,0) înseamnă f la x=1 și y=0.

Afișaj în Mediu Integrat

Sunt folosite trei ferestre:

– o fereastră program unde textul este editat sau selectat

– o fereastră grafică pentru a arăta triangulația curentă

– o fereastră grafică pentru a arăta funcția curentă

Ferestrele grafice pot fi mărite indefinit. Zoom-ul poate fi resetat; atunci afișajul este ajustat astfel încât triangulația se potrivește în fereastră.

Triangulațiile sunt afișate fie direct sau prin micșorarea triunghiurilor; acest lucru vă poate ajuta să verificați dacă nu sunt greșeli.

Frontierele sunt afișate color, în funcție de numerele lor ib.

Funcțiile sunt arătate în trei moduri: linii nivel, hartă color, grafice 1D.

Linii nivel sunt liniile valorilor constante ale funcției.

Harta color este făcută prin colorarea regiunilor dintre două linii de nivel.

Grafic 1 dimensional a unei capturi în domeniu, pe lângă un segment, pe care trebuie să îl definiți cu mouse-ul, apăsând pe primul punct și ținând apăsat mouse-ul până ce atingeți ultimul punct.

Gfem afișează doar câte o funcție, ultima citită sau creată și afișată cu cuvântul cheie grafice ( ).

Programul

Gfem are un limbaj limitat care în general respectă sintaxa limbajului Pascal. Are cuvinte rezervate:

if then else, iter () ( )

dx() dy() id() dnu() laplace(9 convect()

begin end ( de asemenea echivalent pentru „and” )

border (…) buildmesh() onbdy() solve()

plot() savemesh() loadmesh() save () load() exec()

Exemplul 1: Triangulează un cerc și o diagramă f=x*y:

border (1,0,6.28,20) ( x:=cos(t); y:=sin(t) );

buildmesh (200); f=x*y; plot (f);

Exemplul 2: pe cerc / ciclu rezolvă problema Dirichlet -∆u= x*y cu u=0 pe frontiera:

border(1,0,6.28,20) ( x:=cos(t); y:=sin(t) ); buildmesh(200);

onbdy(1) u=0; solve(u) -laplace(u) =x*y; plot(u);

Buildmesh(), border ()

Sintaxa este

border ( ib, tmin, tmax, nbt ) (…x:=f(t); …y:=g(t) …);

… border ( … ) … / * această linie e opțională */

….

buildmesh (nbmax);

Unde f, g sunt funcții generice. Frontiera este dată în formă parametrică. Numele parametrului trebuie să fie „t” iar cele două coordonate trebuie să fie „x” și „y”. dacă parametrul merge de la „tmin” la „tmax”, frontiera trebuie scanată astfel încât să părăsească partea internă a domeniului pe stânga. Frontierele trebuie să fie închise; punctele nbt sunt puse pe aceasta cu valori t = tmin + i * (tmax – tmin) / (nbt-1) / unde i ia valori întregi de la 0 la nbt-1.

Triangulația este creată printr-o metodă Delaunay – Voronoi cu cele mai multe vertice nbmax. Mărimea triunghiurilor este monitorizată de mărimea celei mai apropiate margini de frontieră.

Triangulația poate avea frontiere cu mai multe componente conectate. Fiecare componentă conectată este numerotată de către „ib” întreg.

Frontierele multi – componente sunt definite cu câteva apelări la cuvântul cheie border. Frontiera exterioară trebuie să fie definită invers decât sensul acelor de ceas, dar cu cele interne, fiind găuri, trebuie să fie definite în sensul acelor de ceas datorită faptului că domeniul trebuie întotdeauna să fie pe partea stângă a frontierelor sale orientate.

Obișnuita formulare pascaliana if … then … else poate fi folosită împreună cu formularea compusă: … ( sau început … sfârșit ) (begin/end). Aceasta permite definiții parametrice pe porțiuni ale frontierelor complicate sau poligonale.

Numărul de referință al frontierei ib poate fi suprascris. De exemplu:

border ( 2,0,4,41 )

(

if ( t1 ) then ( x:=t; y:=0 );

if ( (t1) and (t2) ) then ( x:=1; y:=t-1; ib=1);

if ( (t=2) and (t=3) ) then ( x:=3-t; y:=1);

if (t3) then ( x:=0; y:=4-t )

);

buildmesh(400);

aici o parte a pătratului unității are ib=1. Celelalte 3 părți au ib=2.

Variabile geometrice, frontiere interne și regionale

x, y se referă la coordonatele în domeniu

„ib” se referă la numărul de referință frontieră; este zero în interiorul domeniului

„regiune ( region)” se referă la numărul de referință a domeniului care el însuși se bazează pe un număr întreg, ngt, atribuit fiecărui triunghi de către constructorul triangulației.

Este convenabil să existe o frontieră de regiune fie pentru :

– a forța unele margini ale triangulației pe o curbă prescrisă sau pentru

– a separa diferite regiuni astfel ca să li se atribui funcții continue pe porțiuni precum constante dielectrice, modul young, …

Frontierele regionale care se află în interiorul domeniului vor fi numite frontiere interne. Domeniul se află pe ambele părți ale acestora.

Frontierele interne nu sunt închise dar dacă se întâlnesc cu alte frontiere, trebuie să fie la verticele lor. Astfel de frontiere interne pot divide domeniul în câteva sub-domenii.

Fiecare sub-domeniu trebuie să aibă atribuit un număr de regiune. Acest lucru este efectuat de către Gfem, nu de către utilizator. De fiecare dată când border este numit numărul, region este adăugat cu 1. Apoi, ultima margine a fiecărei frontiere interne și externe ( când este invocat cuvântul cheie border ) atribuie numărul region domeniului care se află la stânga sa. Mai exact: De fiecare dată când este evocat cuvântul cheie border, sub-domeniului ( regiunea apropiată definită prin toate frontierele ) de pe partea STÂNGĂ al ultimei margini ale acelei frontiere ( orientată de la următorul la ultimul spre ultimul vertex al acelei frontiere ) îi va fi atribuit un număr de regiune ngt care este egal cu numărul de apariție al cuvântului cheie border înainte în program.

Apoi, de la funcția constantă pe porțiuni ngt, o funcție continuă pe porțiuni region este construită prin atribuirea fiecărui vertex a valori ngt pe triunghiul care este imediat la sud de acesta.

Funcții de construcție

Există patru funcții pre-definite: x, y, ib, region; x, y reprezintă coordonatele ( ca de obicei x este orizontal iar y este vertical ); ib este un întreg care este zero în interiorul și 0 pe frontieră. Pe frontieră este egal cu numărul de referință al fiecărei componente conectate, conform definirii de către border.

Obișnuita formulare if… then… else poate fi folosită cu o restricție importantă la expresia logică care trebuie să returneze o valoare scalar:

if ( expresia logică ) ( formularea; … formulare )

else ( ….. );

„Expresia logică” controlează if – ul dar întoarcerea ei fiind 0 sau 0.

Pot fi folosite auxiliare variabile.

Pentru a minimiza memoria simbolul „:=” spune compilatorului să nu aloce o funcție array acestei variabile. Astfel, pi:=4*antan(1); v=sin(pi*x); generează un array pentru „v” dar nici un array nu este atribuit „pi” – ului.

onbdy și solve

Scopul este de a defini o condiție de frontieră și EDP – urile.

Sintaxa generală este

onbdy (ib) id (u) * exp dnu (u) = g sau onbdy (ib) u = g

unde ib este numărul de frontieră, exp este o expresie generică și g o funcție generică. înseamnă + sau – sau nimic. Termenul id (u) poate fi absent / poate lipsi ca și în –dnu (u) = g; aceasta denotă operatorul de identitate (id(u)=u)

dnu (u) reprezintă conormalul EDP – ului, de exemplu, grad(u).n când operatorul EDP este u-div(.grad(u)).

Sintaxa pentru rezolvare este

solve (u)

id(u)*exp1dx(u)*exp2dy(u)*exp3laplace(u)*exp4=exp5

Aceasta corespunde cu

u exp1 ∂x(u)exp2∂y(u)exp3∆*(exp4∆(u))=exp5.

Observați că laplace(u)* a este o notație care induce în eroare pentru

∂x(a∂xu) + ∂y(a∂yu)

Aceasta rezolvă această EDP prin Metoda de Element Finit de 1 grad despre triunghiuri și o factorizare Gauss. Odată ce matricea este construită și factorizată, solve poate fi invocată din nou prin solve (u-1) …; apoi matricea nu este reconstruită și este efectuată doar o soluție a sistemului linear. Aceasta salvează o grămadă de timpi cpu oricând e posibil. Mai multe matrice pot fi stocate și folosite simultan. În care caz secvența este solve (u,i) …; … solve(u,-i) …; unde „i” este variabil ( nu o funcție array). În orice caz matricele trebuie construite în ordine naturală: i=1 primul, apoi i=2 … apoi ele pot fi folosite în orice ordine. Se poate de asemenea și refolosi o matrice veche, ca și în solve (u,i)…; … solve (u,i) …; solve (u,i) …; Observați că solve (u) este echivalent cu solve (u,1).

Transmitere

Metoda de element finit nu se ocupă cu termenii de convecție propriu – zis, când acestia domină , este necesar upwinding-ul; convect furnizează o unelată pentru upwinding lagrangian pentru că g=convect(f, u, v, dt) este o aproximație a f(x-u(x,y)dt,y-v(x,y)dt). Astfel

∂tf+u∂xf+v∂yf≈[f(x,y,t)-f(x-u(x,y)dt,y-v(x,y)dt,t-dt)]/dt

Aceasta permite schemei să fie implementată ca:

iter(n)

(

Id(f)dt – laplace(f)* m (semn indescifrabil) = h + conect(f, u, v, dt) / dt;

);

Aceasta va rezolva ecuația de convecție – difuzie

∂tf + u∂xf + v∂yf – ∆·(∆f) = h

într-un mod mult mai stabil decât dacă ar fi fost folosiți operatorii dx() și dy().

(amintiți-vă că ∆·(∆f)=∂x(µ∂xf)+∂y(µ∂yf)).

Plot, savemesh, save, saveall, load, loadmesh, exec

În cadrul unui program, cuvântul cheie plot(u) va afișa u.

Instrucția save(’filename’,u) va salva funcția array pe disc.

Instrucția savemesh(’filename’)va salva triangulația pe disc.

În mod similar pentru citirea datelor cu load((’filename’,u) și loadmesh(’filename’). Formatul fișierului poate fi văzut cel mai bine prin deschiderea lor cu un editor text. Pentru o funcție array f, save și load, folosiți:

nb vertices

f[1]

….

f[nb vertices]

Dacă f este o constantă, valoarea sa unică este anexată la sfârșitul fișierului; acesta este de ajutor pentru probleme dependente de timp sau orice problemă cu bucle de ciclare / iterare.

Pentru triangulații loadmesh, savemesh folosiți:

ns nt

q[1].x q[1].y ib[1]

q[ns].x q[ns].y ib[ns]

me[1] [0] me[1][1] me[1][2] ngt[1]

me[nt] [0] me[nt][1] me[nt] [2] ngt[nt]

Se observa că Gfem folosește Fortran standard pentru me[][]și numerotează verticele începând de la numărul 1 în loc de 0 ca și în standard C. însă, în programele C, trebuie să se folosească me[][]-1.

Saveall (varName, ’fileName’) este pentru a salva toate funcțiile care apar în definiția unei EDP. Este destinat pentru aceia care doresc să își scrie propriilor lor solvere. Sintaxa este exact aceeași ca și cea a solve(,), cu excepția faptului că ultimul parametru este numele fișierului. Primul parametru este folosit doar pentru a indica interpretului care este funcția necunoscută. Formatul fișierului este

nbvertices

ng [1] …

ng [i] c[i] g[i] p[i] f[i] b[i] a1[i] a2[i] nu[i]

ng [ns] … … nu [ns]

unde c,g, … corespund EDP

ub+a1∂xu+a2∂yu-∆(·v∆u)=f

sau u=p

Este Dirichlet dacă ng¡0 și întotdeauna Newann / Robin pe rest dar cu coeficient zero dacă utilizatorul nu a definit c și g. ca o regulă, funcțiile nedefinite sunt zero.

Exec (’programe’)

va lansa aplicația program. Aceasta este de ajutor pentru a executa o EDP fortran de exemplu. Pe Macintosh acest lucru face Gfem-ul să se închidă.

Iter

Definește o buclă

Sintaxa:

iter(j) ( … )

unde j se referă la numărul de bucle. Nu sunt permise buclele întrepătrunse.

Condiții de frontieră la colțuri

Colțurile unde condițiile de frontieră se schimbă de la Neumann la Dirichlet sunt chiar mai ambigue datorită faptului că condițiile Dirichlet sunt atribuite verticelor, în timp ce condițiile Neumann trebuie atribuite marginilor de frontieră; Gfem încă nu numără marginile. Pentru condițiile Neumann, Gfem ia la margini înțelesul valorilor la cele două capete de vertice. Dacă un vertex este Dirichlet, atunci valorile unei condiții Neumann nu sunt definite la verticele Dirichlet astfel că e folosit zero.

În cazul unui colț Dirichlet – Neumann este de obicei uzual a se pune un număr Dirichlet la colț. Când acesta este efectuat prin program, trebuie avut în vedere că Gfem verifică ca verticele frontierei să nu fie generate de două ori astfel ca doar prima împrejurare contează.

Derularea & Execuția unei porțiuni de program ( doar mediu integrat )

Comanda ’run ’ va executa întregul program scris în fereastra ’edit’. Același lucru poate fi efectuat prin selectarea întregului text (comanda A) și apăsînd tasta ’enter ’.

O porțiune a programului poate fi executată selectând-o și apăsînd ENTER. Pentru a evita erorile, când un text este selectat în fereastra Edit, toate tastele sunt dezactivate, cu excepția tastei „enter”.

Această caracteristică este de foarte mare ajutor pentru a evita reconstrucția triangulației dacă este o dată deja făcută.

Unele logici trebuiesc aplicate când se folosește această caracteristică altfel vor fi generate mesaje de eroare. Iată regulile:

– textul selectat trebuie să aibă sintaxa corectă a întregului program ( aveți grijă să se potrivească begin și end … )

– Toate funcțiile și variabilele trebuie să fie definite în textul selectat ( și nu în exterior pentru că Gfem nu salvează date, cu excepția triangulației ).

Grafice unidimensionale

Ultima funcție definită de cuvântul cheie plot sau prin deschiderea unui proiect poate fi vizualizată în câteva moduri, dintre care unul fiind o diagramă unidimensională alături de oricare segment definit prin mouse.

Selectarea acestui meniu aduce în față fereastra de triangulație. Gfem așteaptă după input-ul utilizatorului care ar trebui să fie segmentul de linie / liniar pe care funcția trebuie a fi afișată. Astfel, trebuie apăsat mouse – ul la punctul de începere apoi trageți mous – ul și dați-i drumul la punctul final.

Este mai sigur să apăsați în fereastră după inducerea unei redesenări.

Ceea ce se vede este t → f(x(t), y(t)) unde [x(t), y(t)]t este segmentul desenat de către utilizator și f este ultima funcție afișată în fereastra diagramei.

Abscisa este distanța de la începutul punctului segmentului.

Changewait

Acest cuvânt cheie declanșează caracteristica de așteptare – clic care dă posibilitatea utilizatorului să se uite la grafice. În mod normal, când un grafic este afișat, trebuie să se dea clic cu mouse-ul pentru a ieși din execuția programului. changewait modifică această opțiune de a nu da clic sau de a da clic, în funcție de starea anterioară.

Alte meniuri ( doar mediu integrat )

Salvați ca Text se aplică la fereastra frontală.

Salvați ca Postcript este o caracteristică limitată pentru a salva imagini care pot fi incluse intr-un fișier TeX și redimensionat după dorință. Plain Text nu este de foarte mare ajutor pentru poziționarea de imagini într-un text. Astfel postcript-ul generat de către Gfem conține informații pentru un TeX macro care este în fișierul f3d. culorile și numerele se pierd odată cu această caracteristică. Rețineți că este o idee bună să redimensionați fereastra Dvs. la mărimea imaginii Dvs. pentru că va fi inclusă într-un dreptunghi care mărimea ferestrei.

Atribuirea de nb frontieră ( doar Macintosh )

Acest meniu vă permite să modificați numărul de referință frontieră ib.

Trebuie să dați clic pe o margine apoi dați clic pe o altă margine apoi dați clic pe un număr color; toate marginile frontierei, explorate în sensul invers de mers al acelor de ceas ( sau în sensul acelor de ceas dacă este o gaură ) între prima și a doua selecție, va fi atribuit acest număr.

Selecția poate fi modificată cu săgețile stânga, dreapta, sus și jos; stânga și dreapta pentru un capăt, sus și jos pentru celălalt capăt.

Dacă doriți să știți numărul atribuit ultimului sau primului vertex al unei selecții, trebuie să efectuați o redesenare ( comanda Z ). Cu toate că utilizatorul selectează marginile numerele sunt atribuite verticelor, astfel rezultă ambiguitatea. Puteți, de asemenea, să folosiți meniul „Show Numbers” ( arată numerele ).

Pentru afișare, regula este că marginile sunt pictate cu culoarea celui mai mic ib din cele două capete ale verticelor sale.

Meniul „Define Scales” ( doar Macintosh )

Este de ajutor în acest context să reajustați imaginile desenate de mână la mărimea lor fizică. Pătratul cu latura 1 este definit cu mouse-ul prin tragerea lui de sus stânga spre jos dreapta.

Triangulație cu mouse-ul ( doar Macintosh )

MacGfem poate triangula poligoane cu mâna într-un mod rapid și facil:

desenați unul sau mai multe poligoane apropiate dar care nu se intersectează dând clic pe unghiurile poligonului.

– pentru a închide un poligon, dați clic foarte aproape de primul punct sau clic lângă acesta, cu tasta cap apăsată.

În această etapă, puteți adăuga frontiere interne atâta timp cât acestea nu se ating de cele din afară. Aceasta este de ajutor pentru a controla densitatea punctelor.

Puteți, de asemenea, să specificați anumite puncte să devină vertice sau să potriviți domeniul Dvs. prin mișcarea celor date.

Apoi, alegeți o densitate de triangulație din meniu și așteptați.

Automat, domeniul are mărimea ecranului Macintosh. Astfel, x poate varia de la 0 la 640… Este recomandat să schimbați scala folosind meniul potrivit și desenând cu mouse-ul noul pătrat cu latura 1.

Exemple de triangulații

/* cercul ( ring ) unității ( radiusul intern este 0.25) */

twopi := 8 * atan(1);

border(1,0,twopi,60) ( x := cos(t); y := sin(t) );

border(2,0,twopi,20) ( x := 0.25 * cos(-t); y := 0.25 * sin(-t) );

buildmesh(400)

Dreptunghiul [ (0,0), (0,2), (2,1), (0,10) ] */

border (1, 0, 2, 20) ( x:= t; y:= 0 );

border (1, 0, 1, 10) ( x:= 1; y:= t );

border (1, 0, 2, 20) ( x:= 2 – t; y:= 1 );

border (1, 0, 1, 10) ( x:= 0; y:= 1 – t );

/* Un pătrat cu laturi și control bine identificate asupra ib la colțuri */

border (1, 0, 4, 41)

(

if (t1) then ( x:=t; y:=0 );

if ( (t1) and (t2) ) then ( x:=1; y:= t–1; ib:=2 );

if ( (t=2) and (t=3) ) then ( x:= 3–t; y:=1; ib:=3 );

if (t3) then ( x:=0; y:=4-t; ib:=4 );

);

buildmesh(400)

/* Cerc multi-regional */

pi:=4*atan(1)

border(1,0,2,17) ( x:= cos (pi*t); y:= sin (pi*t) );

border(0,-1, 1, 7) ( x:= t; y:= 0 );

border(0, 0, 1, 4) ( x:= 0; y:= t );

buildmesh(300);

Exemple complete

CONDENSOR ELECTROSTATIC

pi:=4*atan(1); /* un cerc cu radiusul 5 centrat la (0,0) */

border(1,0,2*pi,60) begin x:= 5 * cos (t); y:= 5 * sin (t) end;

/* Dreptunghiul pe dreapta */

border(2,0,1,4) begin x:= 1+t; y:= 3 end;

border(2,0,1, 24) begin x:= 2; y:= 3-6*t end;

border(2,0,1, 4) begin x:= 2-t; y:= -3 end;

border(2,0,1, 24) begin x:= 1; y:= -3+6*t end;

/* Dreptunghiul pe dreapta */

border(3,0,1, 4) begin x:= -2+t; y:= 3 end;

border(3,0,1, 24) begin x:= -1; y:= 3-6*t end;

border(3,0,1, 4) begin x:= -1-t; y:= -3 end;

border(3,0,1, 24) begin x:= -2; y:= -3+6*t end;

buildmesh(800);

/* Condiții frontieră și EDP */

onbdy (1) v = 0;

onbdy (2) v = 1;

onbdy (3) v = -1;

solve(v) – laplace (v) = 0;

plot (v);

Transmiterea și radiația căldurii

border(1, 0, 22, 89 )

begin if (t=10) then begin x:= t; y:=0; ib:=3 end;

if (t10) and (t11) then begin x:= 10; y:=t-10; ib:=2 end;

if (t=11) and (t=21) then begin x:= 21-t; y:=1; ib:=4 end;

if (t21) then begin x:= 0; y:= 22 – t; end;

end;

buildmesh(800);

changewait; t0 := 10; tl := 100; te : = 25; b=0.1; c = 5.0e-8;

w = (b + 2*c * (te+546)*(te+273)* (te + 273)) ;

onbdy (1) v=t0; onbdy(2) v = tl; onbdy(3) dnu(v)=0;

onbdy (4) id(v) * w + dnu(v) = te * w;

solve (v,l) -laplace(v) * y =0;

iter (l0)

begin u=v;

w = (b + c * (u + te + 546)*((u+273)*(u+273) + (te+273)*(te+273)));

onbdy(l) v=t0; onbdy(2) v = tl;

onbdy(3) dnu(v)=0; onbdy(4) id(v)*w + dnu(v) = te * w;

solve (v,-1) -laplace(v) * y =0; plot(v);

end

CĂLDURĂ: material non omogen

r0 := 1.0; rl := 2,0;

border (1,0,22,89) begin region :=l;

if (t<10) then begin x:= t; y:=0 ; ib:=3 end;

if ((t>=10) and (t<=ll)) then begin x:=10; y:=rl*(t-10) ;

ib:=2 end;

if ((t>11) and (t<21)) then begin x:=21-t; y:=rl; ib:=4 end;

if (t>=21) then begin x:=0; y:=rl*(22-t) end;

end;

border(0,0,10,41) begin x:= t; y:=r0 end;

buildmesh(800);

t0 = 10; tl = 100; te := 25; kappa =0.01 +max(y-1,0)/(y-1.0001);

onbdy(1) v=t0; onbdy(2) v = tl; onbdy(4) dnu(v)=0.2; onbdy(3)

dnu(v)=0;

solve(v) –laplace(v)*kappa*y +id(v)*kappa*y =0; plot(v);

Flux potențial compresibil

changewait; /* gamma = 1.4, radiusul cercului extern este 5 */

mach1 := l/sqrt(6); machinfty = 0.85*machl;

rhoinfty=sqrt((l-machinftyˆ2)ˆ5);

onbdy(1) dnu(phi) = rhoinfty*machinfty*x/5; onbdy(2) dnu(phi) = 0;

solve(phi) id(phi)*0.0001-laplace(phi) = 0;

ul = dx(phi); u2 = dy(phi); rho=sqrt((l-(ulˆ2+ u2ˆ2)) ˆ5);

plot(phi);

iter(5)

begin

onbdy(l) dnu(phi) =rhoinfty*machinfty*x/5; onbdy(2)

dnu(phi) = 0;

solve(phi) id(phi)*0,0001-laplace(phi)*rhop=0;

ul = dx(phi); u2 = dy(phi); rho=sqrt((1-(ulˆ2+ u2ˆ2))ˆ5)

rhop = convect(rho,ul,u2,0.1); plot(rho)

end;

plot(sqrt((ulˆ2+u2ˆ2))/machl);

Ecuații Navier – Stokes

/* redus / slab dar mai bune decât nici un algoritm */ changewait;

border (1,0,1,6) begin x: =0; y:=l-t end;

border(2,0,1,15) begin X:=2*t; y:=0 end;

border (2,0,1,10)begin x:=2; y: = -t end;

border(2,0,1,20) begin x:=2+3*t; y.=-l end;

border(2,0,1,35)begin x:=5+15*t; y:=-l end;
border(3,0,1,10) begin X:=20;y:=-l+2*t end;

border(4,0,l,35) begin x:=5+l5*(1-t); y:=l end;

border(4,0,1,40) begin x:=5*(1-t); y:=1 end;
buildmesh(900);

nu = 0.002; dt := 0.4;

/* presiunea inițială */

onbdy(1)dnu(p) =-2*nu; onbdy(2) dnu(p) = 0;

onbdy(3) p=0; onbdy(4) dnu(p) = 0;

solve(p,1) – laplace(p)= 0;

/* viteza inițială orizontală */

onbdy(l) u = y*(l-y); onbdy(2) u = 0;

onbdy(3) dnu(u) = 0; onbdy(4) u = 0; solve(u,2)id(u)/dt-laplace(u)*nu=-min(y*y-y,0)/dt;

/* viteza inițială verticală */

onbdy(l)v = 0; onbdy(2) v = 0;

onbdy(3) v = 0; onbdy(4) v = 0;

solve(v,3) id(v)/dt-laplace(v)*nu =0; un = u, vn = v; iter (80)

begin f = convect(un,u,v,dt); g = convect(vn,u,v,dt);

/* VITEZĂ orizontală */

onbdy(l) u = y*(l-y); onbdy(2) u = 0;

onbdy(3)dnu(u)=0; onbdy(4) u = 0; ,

solve(u,~2) id(u)/dt-laplace(u)*nu = f/dt -dx(p);

plot(u);

/* VITEZĂ verticală */

onbdy(1) v = 0; onbdy(2) v = 0; onbdy(3) v = 0; onbdy(4) v = 0;

solve(v,-3) id(v)/dt-laplace(v)*nu = g/dt -dy(p);

/* presiune */

onbdy(l)dnu(p) =-2*nu; onbdy(2) dnu(p) = 0;

onbdy(3) p=0; onbdy(4) dnu(p) = 0;

solve(p,-l) -laplace(p)= -(dx(£) + dy(g))/dt;

un = u; vn = v;

end ;

save('u.dta',u); save( 'v.dta',v); save('p.dta',p); plot(u);

Bi-NACA cu flux de ridicare potențial

twopi := 8 * atan(l);

border(1,0,twopi,40)

( x:=5*cos(t); y:=5*sin(t) );

border(2,0,2,71)

(

if (t< = l) then

(

x : = t ;

y := 0.17735*sqrt(t)-0.075597*t- 0.212836*(tˆ2)+ 0.17363*(tˆ3)-0,06254*(tˆ4);

) else

(

x := 2-t;

y:= – (0.17735* sqrt (2-t) -0.075597* (2-t)-0.212836*((2-t) ˆ2)+0.17363*((2-t) ˆ3) -0.06254* ( (2-t) ˆ4))

);

);

border(3,0,2,71)

(

if(t<=l) then

(

x := t+0.5;

y := 0.5+0.17735*sqrt(t)-0.075597*t

– 0.212836* (tˆ2)+0.17363*(tˆ3)-0.06254*(tˆ4); )

) else

(

x := 2-t+0.5;

y:= 0.5-(0.17735*sqrt(2-t)-0.075597* (2-t)-0.212836

*((2-t)ˆ2)+0.17363*((2-t)ˆ3)-0.06254*((2-t)ˆ4))

);

);

buildmesh(2800) ;

onbdy(l) psi0 = y-0.1*x; onbdy(2) psi0 = 0;

onbdy(3) psi0 = 0;

solve(psi0) -laplace(psi0) = 0; plot(psi0);

onbdy(l) psil = 0; onbdy(2) psil = 1;

onbdy(3) psil = 0;

solve(psil,-1) -laplace(psil)= 0; plot(psil);

onbdy(l) psi2 = 0; onbdy(2) psi2 = 0;

onbdy(3) psi2 = 1;

solve(psi2,-1) -laplace(psi2)= 0; plot(psi2);

beta := psi0 (0.99,0.01)+psi0(0.99,-0.01);

beta := -beta / (psil(0 .99,0.01)+ psil(0.99,-0.01)-2) ;

beta2 := psi0(1.49,0.51)+psi0(1.49,0.49);

beta2 := -beta2 / (psi2(1.49,0.49)+ psi2(1.49,0.51)-2) ;

psi = beta*psil +beta2*psi2+psi0; plot(psi);

cp = -dx(psi) ˆ2 – dy(psi) ˆ2; plot(cp);

/* obiectul este de a găsi o condiție de frontieră eficientă pentru a evita acest tip de calculație care a fost efectuată aici pentru a verifica rezultatul. */

(

pi := 4*atan(l);

teta := pi/4;

re := 10;

ri := 3;

c := 9.5;

rr:= 0.30;

border(1,0,15,300)

( t3:=12-t;

if (t3<=l) then

(

t1 := teta + (3+t3)*(pi/2-teta)/4;

y := re * cos (t1) ;

x := re * sin(t1);

) else

if ((1<t3) and (t3<2)) then

(

t1 := (t3-l)*(re-c-rr);

y := 0;

x := re – t1;

) else

if ((2< = t3)and(t3< = 3)) then

(

t1 := (t3-2)*pi;

y : = rr * sin (t1) ;

x := c + rr * cos (tl) ;

) else

if ((3<t3) and (t3<6)) then

(

t1 := (t3-3) * (c-rr-ri)/3;

y := 0;

x : = c – rr – tl;

) else

if ((6<=t3) and (t3<=7)) then

(

tl := pi/2+(t3-6)*(teta-pi/2) ;

ib := 2;

y := ri * COS (tl) ;

x := ri * sin(tl);

) else

if (7<t3) and (t3<10) then

(

tl := ri + (t3- 7) * (c – ri-rr) /3;

ib := 3;

y := tl * cos(teta) ;

x := tl * sin(teta);

)

else

if (t3>=10) and (t3<=11) then

(

tl :=-3* pi/4+(t3-10)*pi;

y := c*sin(pi/4) + rr * sin(tl);

x := c*cos(pi/4) + rr * cos(tl);

) else

if (11<t3) and (t3< = 12) then

(

tl := c+rr+(t3-11)*(re-c-rr);

ib := 3;

y := tl*sin(pi/4);

x := tl*sin(pi/4);

)

);

border(1, 0, 2, 30)

(

tl:=-t*pi;

y :=c*sin(pi/8) + rr * sin(tl);

x := c*cos[pi/8) + rr * cos(t1);

);

border(l, 0, 2, 30)

(

tl :=- t*pi;

y :=c*sin(pi/16) + rr * sin(tl);

x := c*cos(pi/16) + rr * cos(tl);

);

border(1, 0, 2, 30)

(

tl := – t*pi;

y :=c*sin (3*pi/16) + rr * sin(tl);

x := c*cos(3*pi/16) + rr * cos(tl);

);

border (l, 0, 2, 30)

(

tl := -t*pi;

y := c*sin (pi/8+pi./32) + rr * sin(tl);

x := c*cos (pi/S+pi/32) + rr * cos(tl);

);

border (1, 0, 2, 30)

(

tl := – t*pi;

y := c*sin (pi/4-pi/32) + rr * sin(tl);

x := c*cos (pi/4-pi/32) + rr * cos(tl);

);

border (l, 0, 2, 30)

(

t1 :=- t*pi;

y := c*sin (pi/32) + rr * sin(tl);

x := c*cos (pi/32) + rr * cos(tl);

);

border (1, 0, 2, 30)

(

tl := – t*pi:

y := c*sin (pi/8-pi/32) + rr * sin(tl);

X := c*cos (pi/8-pi/32) + rr * cos(tl);

);

buildmesh(3400);

savemesh('wallBig.msh');

onbdy(l) v a 0;

onbdy(2) v = 0;

onbdy(3) dnu (v) = 0;

solve(v) -laplace(v)=x*y/sqrt(x^2ty^2);

plot(v);

save('wlb.dta' ,v) ;

a:=l.4;

r:=8;

11:= a*r/ 110-r);

nu := ((r*(3 + 11)-3*(2 + ll)) / ((11-2)/(r^4) – (2 + 11)/8l))/10;

mu := 3/10 – nu/81;

u = 2*x*y*(mu + nu/((x^2+y^2)^2) – sqrt(x^2+y^2)/10);

load('wlb.dta', v) ;

plot(v-u);

II.4 Studii de caz

Abstract

Vom ilustra cu câteva exemple, că multe probleme de fizică, matematică și altele pot fi formulate cu Ecuații diferențiale Parțiale ( EDP ) și rezolvate cu interpretul Gfem-ului.

II.4.1 ELECTROSTATICA

Formularea problemei

Fie {Ci}1 … N iar N conductorii într-o îngrădire, C0 , fiecare întreținută la un potențial electrostatic φi. Să presupunem că îngrădirea este un conductor de potențial 0. dacă dorim să știm potențialul electrostatic φ(x) la un punct x, trebuie să rezolvăm:

∆φ = 0 în Ω φ|Γ=g

unde Ω este interiorul îngrădirii minus conductorii Ci și Γ este frontiera lui Ω, adică, îngrăditura plus suprafața obiectelor Ci. Aici g este orice funcție a lui x care este egală cu φ pe Ci și 0 pe C0. astfel, a doua ecuație este o notare scurtă pentru

φ = φi pe Ci i=1 … N, φ = 0 pe C0

amintim din nou definiția operatorului Laplace ∆:

Aceasta este numită o problemă dirichlet pentru ecuația Laplace. Poate fi arătat ( Dautray – Lions [1989], Godlewski – Raviart [1992], Strang [1989]) că acesta are o soluție unică.

Astfel de probleme reies în proiectarea condensatorilor, printre altele.

Fig. 3

Model de condensator realizat din doi electrozi

dreptunghiulari paraleli într-o îngrăditură circulară.

Observați simetria din această problemă, în ceea ce privește axele verticale și orizontale trecând prin centrul cercului.

Triangulația

Cercul

Domeniul este creat din trei frontiere: un cerc și două dreptunghiuri. Ecuația cercului în formă parametrică este:

x(t) = R cos t y(t) = r sin t, t ε ]0,2π]

În limbajul Gfem aceasta devine:

Pi:= 4*atan(1);

border (1, 0, 2*pi,60)

begin

x:=5*cos(t);

y:=5*sin(t);

end;

informația adițională este aceea că

– cercului i s-a dat un număr de referință de frontieră egal cu 1.

– radiusul său este R = 5.

– parametrul t îl scanează în sens opus acelor de ceasornic astfel că domeniul din interiorul cercului va fi triangulat, nu în exterior / afară ( dar exteriorul fiind independent Gfem ar fi detectat o eroare )

– există 60 de puncte pe cerc care vor deveni vertice ale triangulației. De fapt sunt generate 61 de puncte dar Gfem recunoaște că primul și ultimul punct sunt aceleași. Astfel trebuie introduse frontierele închise.

Dreptunghiurile

Poligoanele pot fi introduse în două moduri:

– ca un set de segmente drepte fiecare descris de un bloc de frontieră; de exemplu dreptunghiul [(1,-3); (2,-3); (2,3); (1,3)] este

border(2,0,1,4) begin x: = l + t; y:=3 end ;

border (2,0, 1,24) begin x:=2; y:=3-6*t end ;

border(2,0,1,4) begin x:=2-t; y:=-3 end ;

border(2,0,1,24) begin x:=l; y:=-3+6*t end

informația suplimentară fiind că este frontiera numărul 2 și fiecare parte are respectiv 4, 24, 4, 24 vertice. Din moment ce este scanat în sensul de mers al acelor de ceasornic, triangulația va fi în afară.

– A doua metodă este de a descrie poligoanele doar printr-un singur bloc frontieră.

border ( … ) begin …

if ( (ti<t) and (t <= tiplusone )) then

x:= … y:= …

else if

… end;

Haideți să alegem prima metodă. Trebuie avut grijă ca fiecare capăt vertex al unui segment să coincidă cu începutul vertexului următorului segment.

Programul final este :

/* un cerc cu raza 5 centrat la (0,0) */

border(l,Q,2*pi,60)

begin

x := 5 * cos (t);

y:= 5 * sin ( t );

end;

/* Dreptunghiul de pe dreapta */

border(2,0,1,4) begin x;=l+t; y:=3 end ;

border(2,0,1,24) begin x:=2; y:=3-6*t end ;

border(2,0,1,4) begin x:=2-t; y:=-3 end ;

border(2,0,1,24) begin x:=l; y:=-3+6*t end ;

/* dreptunghiul de pe stânga */

border(3,0,1,4) begin x:=-2+t; y:=3 end ;

border(3,0,1,24) begin X:=-l; y:=3-6*t end ;

border (3,0,1,4) begin x:=-l-t; y;=-3 end ;

border(3,0,1,24) begin X:=-2; y:=-3+6*t end

buildmesh(800);

/* condiții de frontieră și EDP */

onbdy(l) v = 0;

onbdy(2) v = 1;

onbdy(3) v = -1;

solve(v) -laplace(v) =0;

plot(v);

Mesh-ul are cel mult 800 de puncte și toate frontierele au condiții Dirichlet. Aici punem 1 volt pe un electrod, -1 volt pe celălalt și 0 pe îngrăditura circulară.

Fig.4

Mesh 800 puncte

O alternativă pentru poligoane

O descriere parametrică a unui poligon de vertice [(a, a)] este

t ε ]ti, ti+1]

Pentru a alege ti-ul trebuie luat în calcul distribuția punctelor pe frontieră. într-adevăr punctele vor fi uniform distribuite între tmin și tmax, astfle că dacă doriți de două ori atâtea puncte pe segment aj+1 – aj decât pe ai+1 -ai atunci trebuie să avem

tj+1-tj=2(ti+1-ti)

în final, din moment ce dreptunghiurile sunt găuri în domeniu, acestea trebuie scanate în direcția acelor de ceas.

De exemplu, un dreptunghi cu vertice [1,-3], [2,-3], [2,3], [1,3] și 56 puncte ar fi

border(2,0,14,57)

begin

if(t<=l) then begin x := 2 – t; y :=-3 end; if((l<t)and(t<=7)) then begin x:=l; y:=-4+t end; if ((7<t)and(t<=8) ) then begin x:= t-6; y:=3 end else if(8<t) then begin x:=2; y:=ll- t end

end ;

Alegerea 57 este justificată ca 14*4+1 dând 4 puncte per unitate, crescând cu t astfel că colțurile devin vertice.

Cel de-al doilea dreptunghi este obținut din primul prin traducerea lui -2 înx.

Programul final

pi := 4*atan (l);

border(1,0,2*pi,60)

begin

x:=5*cos(t);

y:=5*sin(t);

end;

border(2,0,14,57)

begin

if (t<=l) then begin x: = 2 – t; y:= -3 end

if((l<t)and(t<=7))then begin x:=l;y:=-4+t end;

if((7<t)and(t<=8)) then begin x:=t-6;y:=3 end

else if (8<t) then begin x:=2; y:=11-t end

end;

end;

border(3.0,14,57)

begin

if (t<=l) then begin x := -1-t; y:= -3 end; if ((l<t) and(t<=7))then begin x:= -2; y:=-4+t end; if ((7<t) and (t<=8))then begin x:= t-9; y:= 3 end else if(8<t) then begin x:= -1; y:=ll-t end

end

buildmesh (800)

onbdy(l) v = 0;

onbdy(2) v = 1;

onbdy(3) v = -1;

solve(v) -laplace(v) =0;

plot(v)

II.4.2 MEMBRANELE

Expunerea problemei

O membrană elastică Ω este lipită de un suport drept și rigid Γ dar greutățile apasă pe suprafața sa ( de exemplu, un lichid ) astfel că pe fiecare parte de suprafață elementară dx=dx1dx2 există o greutate f(x)dx. Atunci, așezarea verticală φ(x) a membranei este soluția pentru

– ∆φ = f în Ω φ|Γ = 0

Fig.5

a.Membrana b.Diagrama

Figura 5: a.- este arătată o membrană din perspectivă, fixată pe o ramă eliptică ( linia îngroșată ) și deviată de o forță de suprafață f.

b.- o diagramă unidimensională a soluției pe axa orizontală

Membrana eliptică

O elipsă a radiusului orizontal a=2 și radiusul vertical 1 și o forță de suprafață egală cu 1 ar da

pi : = 4 * atan(1);

a : = 2;

border(1.0,2*pi,40)

begin

x := a * cos(t);

y : = sin(t);

end ;

buildmesh(800) ;

onbdy(1) v = 0;

solve(v) -laplace(v) = 1;

plot(v);

Analize de eroare

Această problemă are o soluție analitică, astfel că putem studia eroarea ca pe o funcție a numărului verticelor.

Când a=1 soluția analitică este:

Astfel că pentru diagramă adăugăm pur și simplu,

plot(v-(1-x*x-y*y)/4);

II.4.3 TRANSMITEREA CĂLDURII

3.1 Grindă rectangulară

O grindă lungă în formă de cărămidă este încălzită la ambele capete astfel că un capăt este la temperatura de t0 iar celălalt la t1> t0, dar pierde căldură prin părțile sale pentru că aerul este la temperatura de te<ti. Presupunând că grinda este suficient de lungă și subțire, problema este formată în două dimensiuni și pierderea de căldură de către părți este formată de către disipare α proporțională cu t0+t1-2te. când grinda Ω este omogenă, ecuațiile pentru temperatura u sunt:

-∆u + αu = 0 în Ω

u|Γ1 = t0, u|Γ2 = t1, |Γ3 = 0

unde Γ = ∂Ω, Γ1 este un capăt, Γ2 celălalt capăt iar Γ3 = ∂Ω – Γ1 – Γ2 cele două părți laterale.

Pentru o grindă de lungime 10 programul este

border(1,0,22,89)

begin

if (t< = 10) then

begin x:= t; y:=0 ; ib:=3 end;

if((t>10) and (t<11))then

begin x:=10; y:=t-10; ib:=2 end;

if ((t>=11) and (t <=21)) then

begin x:=21-t; y:=l; ib:=3 end;

if (t>21)then

begin x:=0; y:=22-t end;

end;

buildmesh(800);

t0 := 10; tl := 100; te := 25;

alpha =0.1 * ( (t0 + tl)/2 – te);

onbdy(l) v=t0;

onbdy(21 v = tl;

onbdy(3) dnu(v)=0;

solve(v) -laplace(v) + id( v ) * alpha = 0;

plot(v);

Grindă axisimetrică

Dacă grinda are o captură circulară atunci problema poate fi stabilită în trei dimensiuni și încă rezolvată în două, profitând de axisimetrie.

În coordonatele cilindrice operatorul Laplace este

Din moment ce nu există nici o dependență θ obținem următoarea problemă

-∂r(r∂ru) – ∂z(r∂zu) + rαu = 0 în Ω

u|Γ1 = t0, u|Γ2 = t1, |Γ3 = 0

Observați șiretlicul prin care ecuația a fost multiplicată cu r astfel încât să se obțină același coeficient pentru primii doi termeni; aceștia sunt de asemenea egali cu -. Partea de triangulație fiind aceeași, acum programul pentru partea de EPD este

t0 := 10; tl := 100; te := 25;

alpha =0.1 * y * ( (t0 + tl)/2 – te);

onbdy(l) v=t0;

onbdy(2) v = tl;

onbdy(3) dnu(v)=0;

solve(v) -laplace(v)*y *id( v ) * alpha =0;

plot(v);

3.3 Convecție și iradiere

Pierderea de căldură prin convecție și iradiere a fost considerata până acum ca o disipare α. Dar o modalitate mai bună de a modela răcirea grinzii prin aerul exterior este folosirea unei pierderi de căldură ale părților proporțională cu diferența dintre temperatura peretelui și temperatura exterioară. Aceasta dă :

Pentru iradiere, legea obiectelor negre dă o pierdere de căldură proporțională cu temperatura absolută la puterea a patra:

Astfel, condiția generală pentru convecție și iradiere este

Această problemă este non-lineară astfel că trebuie rezolvată iterativ. Dacă m denotă indexul iterativ, atunci ne putem aștepta să rezolvăm problema prin setarea pe laterală în contact cu aerul condiția de frontieră.

Această condiție vine de la identitatea

a4 – b4 = ( a – b ) ( a + b ) ( a+b2 )

border(1,0,22,89)

begin

if(t<=10)then

begin x:= t; y:=0 ; ib:=3 end;

if((t>10)and(t<ll))then

begin x:=10; y:= t-10; ib:=2end;

if ((t>=ll)and(t<=21))then

begin x:=21-t; y:=l; ib:=4 end;

if(t>21)then

begin x:=0; y:=22-t end;

end;

buildmesh(800);

Changewait;

t0:=10; t1:=100; te:=25;

b=0.1; c=5.0e-8

w=(b+2*c*(te+546)*(te+273)*(te+273));

onbdy(1) v=t0;

onbdy(2) v=t1;

onbdy{3) dnu(v)=0;

onbdy(4) id(v) * w + dnu(v) = te * w;

solve(v,l) -laplace(v) * y = 0;

iter(10)

begin

u=v;

w=(b+c*(u+te+546)*(u+273)*(u+273)+(te+273)*(te+273)));

onbdy(l) v=t0;

onbdy(2) v = tl;

onbdy(3) dnu(v)=0;

onbdy(4) id(v)*w + dnu(v)= te * w;

solve(v,-l) -laplace(v) * y =0;

plot(v);

end

Fig.6 Pierderea de caldura

3.4 Mediu non omogen

Luați în considerare o grindă executată din două materiale diferite. Materialul i ocupă Ωi și are iradierea termică ki. domeniul pe care temperatura este studiată este Ω făcut prin unirea Ω1 cu Ω2.

Când se ajunge la un echilibru și pentru o grindă axisimetrică, temperatura u satisface EDP

în Ω

Această ecuație conține condiția de discontinuitate asupra fluxurilor de căldură de-a lungul ∑, frontiera obișnuită între Ω1 și Ω2.

Totodată, κ(x) este egal cu κi când x Ωi. Amintim definiția operatorilor div și grad:

Să presupunem că grinda este făcută din două grinzi coaxiale, una cu o secțiune transversală circulară de raza r0 și cealaltă cu un inel precum o secțiune transversală cu raza r0 și r1.

Pentru a triangula un astfel de domeniu, vom folosi o frontieră artificială la interfața între cele două tipuri de material. Procedând în acest mod, trebuie să avem grijă ca punctele sale de capete să coincidă cu verticele celeilalte frontiere. De asemenea, orientarea trebuie să fie de așa natură, ca fiecare subdomeniu să aibă una din ultimul capăt al fiecărei frontiere pe partea sa stângă. Într-un subdomeniu setăm regiunea = 1 astfel că această variabilă ar putea fi folosită mai târziu pentru a defini coeficienții discontinuu. Fiecărui triunghi îi este atribuit un întreg pentru valoarea regiunii.

r0:=1.0; r1:=2.0;

border(l,0,22,89)

begin

region :=1;

if(t<10)then begin x:= t; y:=0 ; ib:=3 end;

if((t>=10)and(t<=ll))then

begin x:=10; y:=rl*(t-10); ib:=2 end;

if((t>ll)and(t<21))then begin x:=21-t;y:=rl; ib:=4 end;

if(t> = 21)then begin x:=0; y:=rl*(22-t) end;

end;

border(0,0,10,41)

begin x:= t; y:=r0 end;

buildmesh(800);

t0 = 10; tl= 100; te := 25;

kappa =0.01 + max(y-1,0)/(y-1.0001);

kappa = y * kappa;

onbdy(l) v=t0;

onbdy{2) v = tl;

onbdy(4) dnu{v)=0.2;

onbdy(3) dnu(v)=0;

solve(v) -laplace(v)*kappa+id(v) =0;

plot(v);

am luat o difuzie termică egală cu

κ= 0.01 în Ω1, κ = 1.01 în Ω2

și o disipare α=κ. Apoi pentru a face schimbarea materialului mult mai influentă am luat o pierdere de flux de căldură pe o parte egală cu 0,2.

Fig.7 Transferul de caldura

II.4.4 ACUSTICA

Variațiile de presiune în aer în repaus sunt guvernate de ecuația undelor, dar când variațiile sunt periodice în timp cu o singură frecvență κ, o ecuație pentru amplitudine v(x) a presiunii sunetului este suficientă pentru a descrie variațiile:

k2v + c2∆v = f în Ω

unde c este dispersarea sunetului, f sursele de sunet.

Pentru a face un poligon cu dungi triunghiulare, care se aseamănă cu forma unei săli de concerte, trebuie să descriem fiecare din marginile sale printr-o frontieră. Aceasta este mult mai ușor cu condiția ca punctul capăt al fiecăruia să întâlnească punctul de început al următorului.

border(1,0,1,20) begin x:=5; y:=1+2*t end;

border(1,0,1,20) begin x:=5-2*t; y:=3 end;

border(1,0,1,20) begin x:=3-2*t; y:=3-2*t end;

border(1,0,1,10) begin x:=1-t; y:=1 end;

border(1 0,1,10) begin x:=0; y:=1-t end;

border(1,0,1,10) begin x:=t; y:=0 end;

border(1,0,1,20) begin x:=1+4*t; y:=t end;

buildmesh(1000);

f=exp(-10*((x-0.5)*(x-0.5)+(y-0.5)*(y-0.5)));

onbdy(1) dnu(v)=0;

solve(v) laplace(v)+id(v)*0.8=f

plot(v);

în zona de etapă, atribuim o sursă descompusă exponențială .

Lungimea undei și viteza sunetului au fost în așa mod alese pentru a se observa unele rezonanțe.

Fig.8 Variatia de presiune

II.4.5 CURGERE FLUIDE

Curgere irotationala fara frecare

Curgeri incompresibile

Când sunt neglijate efectele viscozității și nu există vârtejuri într-o curgere plană paralelă, un exemplu bun de viteză a fluidului u(x) la punctul x= (x1,x2) este dat de:

unde ψ satisface ∆ψ = 0

Această ecuație se traduce prin faptul că fluidul este incompresibil (·u=0) precum apa sau aerul la viteză moderată, și că nu există vârtejuri (xu=0). Fluidul nu penetrează pereți solizi astfel u trebuie să fie tangent cu pereții solizi. Tradus în funcția curgerii ψ aceasta este:

ψ este constant față de pereți / pe pereți

pentru o țeavă pentru care ratele de aflux și de eflux sunt cunoscute, problema ar fi

-∆ψ = 0 în Ω, ψ|Γ= g

cu ψ linear pe frontierele de aflux și eflux și ψ constant la pereți.

Luați în considerare o duză care este simetrică în ceea ce privește cele două axe x,y și care trece prin punctele ( 0,1 ) și ( a, 1 + b ) cu înclinare orizontală. ecuația sa ar fi

-a < x< a

cu o viteză unitară la frontiera de intrare problema este rezolvată prin

a:= 6; b:= 1;

border(1,0,1,8)beginx:=-a;y:=1+b-2*(1+b)*t end;

border(2,0,1,16) begin x:=-a+2*a*t;

y:=-1-b*(x/a)*(x/a)*(3-2*abs(x)/a)

end;

border(3,0,1,8)beginx:=a;y:=-1-b+2*(1+b)*t end;

border(4,0,1,16) begin x:=a-2*a*t;

y:=1+b*(x/a)*(x/a*(3-2*abs(x)/a)

end;

buildmesh(800);

onbdy(1) psi=y;

onbdy(2) psi=-1 –b;

onbdy(3) dnu(psi)=0;

onbdy(4) psi=1+b;

solve(psi) –laplace(psi)=0;

plot(psi);

Fig.9 Curgere irotationala

5.2 PROFILE LIFTING

în orice caz, în mai multe situații complexe valoarea constantelor pentru ψ se poate să nu fie cunoscută și se poate să trebuiască să fie calculată ca și în cazurile profilurilor cu ridicare. haideți să studiem caracteristicile aerodinamice a unui profil de aripă S în într-un tunel aerodinamic; doar parte din tunelul aerodinamic este parte a unui domeniu computațional și o frontieră care este departe de aripă dar încă în cadrul tunelului aerodinamic este numită Γ∞. În acest caz, pentru curgeri incompresibile, ca și mai înainte,

∆ψ = 0 în Ω

ψ|s = c

ψ|Γ∞ = u∞2x – u∞1y

unde Ω este regiunea de spațiu ocupat de fluid, u∞ este viteza ( dată ) a fluidului departe de profil; c este o constantă care trebuie să fie în așa fel încât ∂nψ este continuu la marginea bordului de scurgere P din S ( valoarea părții de sus a profilului se potrivește cu valoarea de pe partea de jos a marginii bordului de scurgere. atunci ( vezi Landau – Lifschtiz [1962]) ridicătura aripii este proporțională cu c.

pentru a găsi c-ul, trebuie să se folosească suprapuneri. soluția ψc este lineară în c:

ψc = ψ0 + c ( ψ1 – ψ0 )

unde ψ0 este soluția EDP cu c=0 și ψ1 cu c = 1. Astfel că odată ce aceste două sunt calculate c este determinat prin scrierea continuității lui ∂ψ / ∂n la marginea bordului de scurgere.

Exemplu de aplicație: ridicarea unei aripi pentru un planor sau un avion la decolare.

Ecuația suprafeței de sus a aripii NACA0012 este

y:=0.17735√x – 0.075597x – 0.212836×2 + 0.17363×3 – 0.06254×4

Mai întâi construim o frontieră circulară externă pentru a aproxima câmpul depărtat. Apoi, ambele părți ale aripii NACA sunt discretizate ( ib = 2 ). în final, pentru un unghi α astfel că tanα = 0.1, problema rezolvată fiind

-∆ψ = 0 în Ω, ψ|Γ1= y – 0.1x, ψ|Γ2= c

găsim c-ul calculând

-∆ψ0 = 0 în Ω, ψ0|Γ1= y – 0.1x, ψ0|Γ2= 0

-∆ψ1 = 0 în Ω, ψ1|Γ1= y – 0.1x, ψ1|Γ2= 1

Soluția generală fiind ψ = ψ0 + cψ1 găsim c-ul scriind că ∂nψ nu va sări la marginea P = (1,0). Acum, ∂nψ ≈ – (ψ(P+) – ψ(P))/δ unde P+ este chiar deasupra P în direcția normală la distanța δ. Astfel, săritura lui ∂nψ este ψ0|P++ ψ0|P-+c(ψ1|P+ + ψ1|P- ) – 2 ) împărțit la δ datorită schimbării normale a semnelor sale de sus și de jos de P. Astfel

care va fi programat astfel

twopi:=8*atan(1)

border(1,0,twopi,40) begin x:=5*cos(t); y: =5*sin(t) end ;

border(2,0,2,71)

begin

if (t<=l) then

begin

x := t;

y := 0.17735*sqrt(t)-0.075597*t

– 0.212 836*(t*t)+0,17363*(t*t*t)-0.06254*(t*t*t*t);

end

else begin

X := 2-t;

y:= -(0.17735*sqrt(2-t)-O.075597

*(2-t)-0.212836*((2-t)*(2-t))

+0.17363*((2-t)*(2-t)*(2-t))

-0.0 6254*(2-t)*(2-t)*(2-t)*(2-t))

end ;

end ;

buildmesh(1800);

onbdy(l) psi0 = y-0.1*x;

onbdy(2) psi0 = 0;

solve(psi0) -laplace(psi0) = 0;

plot(psi);

onbdy(l) psil = 0;

onbdy(2) psil = 1;

solve(psil,-l) – laplace(psil)= 0;

plot(psil);

beta := psi0 (0.99,0.01) +psi0(0.99,-0.01) ;

beta:= -beta/(psil(0.99,0.01)*+psil(0.99,-0.01)-2);

psi = beta*psil *psi0;

plot(psi);

cp = dx(psi)ˆ2 – dy(psi)ˆ2 ;

plot(cp);

norma gradientului ψ este afișată pentru că este presiunea.

Observați, din nou, folosirea matricei factorizate când rezolvați EDP a doua oară.

Fig.10 Profil aerodinamic aripa

5.3 CURGERE IROTATIONALA COMPRESIBILA

Curgerile adiabatice compresibile au o densitate ρ dată de

ρ = ( 1 – |u|2)

unde, ca și mai înainte, viteza u este aceea că

ρu= (∂2ψ, -∂ψ1)T, ∂2u1 – ∂1u2 = 0

Pe scurt,

cu condiții de frontieră similare ca și în 5.2.

Atâta timp cât curgerea rămâne subsonică peste tot , aceasta poate fi rezolvată prin metoda iterativă a punctului fixat / fix, precum

Partea de triangulație a programului fiind aceeași ca în paragraful 5.2., restul programului este

changewait;

/* gamma = 1.4 */

machl := 1/sqrt(6);

machinfty = 0.6*machl;

onbdy(l) psi = machinfty*y;

onbdy(2) psi = 0;

solve(psi) -laplace(psi) = 0;

u2 = -dx(psi); ul = dy(psi);

rho=sqrt((l-(ulˆ2+ u2ˆ2)) ˆ5);

plot(psi);

iter(5)

begin

onbdy(l) psi = machinfty*y;

onbdy(2) psi = 0;

solve(psi) –laplace(psi)/rho=0;

u2 = -dx(psi)/rho; ul = dy(psi)/rho,

rho=sqrt((l-(ulˆ2+ u2ˆ2)) ˆ5);

plot(rho)

end;

plot(sqrt((ulˆ2+u2ˆ21)/machl);

Acest algoritm nu va funcționa pentru numere mari ale lui Mach, de aceea este mai bine să se lucreze cu o funcție potențială decât cu o funcție de curent. Aceeași problemă este scris ca

ρ=(1-u2)2.5, u= (∂1ø, ∂2ø )

Acum, condițiile de frontieră sunt

O schemă iterativă ca și

este stabilă atâta timp cât numărul Mach rămâne mai mic de 1.

5.4 CURGEREA TRANS-SONICA

dacă este introdus vânt ascendent, atunci numărul mach poate deveni supersonic în unele locuri. Vântul ascendent poate fi obținut cu funcția convect care face

ρ_(x)=ρ (x-u(x)dt)

Acum programul este

changewait;

/*gamma = 1.4,radiusul cercului extern e 5 */

machl := l/sqrt(6);

machinfty = 0.85*machl;

rhoinfty=sqrt((1-machinftyˆ2) ˆ5);

onbdy(1) dnu(phi) = rhoinfty*machinfty*x/5;

onbdy(2) dnu(phi) = 0;

solve(phi) id(phi)*0.0001-laplace(phi) = 0;

ul = dx(phi); u2 = dy(phi);

rho=sqrt((1-(ulˆ2+ u2ˆ2)) ˆ5);

plot (phi) ;

iter(5) begin

onbdy(1) dnu(phi) =rhoinfty*machinfty*x/5;

onbdy(2) dnu(phi) = 0;

solve(phi) id(phi)*0.0001-laplace(phi)*rhop=0;

u1=dx(phi); u2=dy(phi);

rho=sqrt ((1-(u1ˆ2+u2ˆ2)) ˆ5;

rhop=convect(rho, u1, u2, 0.1);

plot (rho)

end;

plot(sqrt((u1ˆ2+u2ˆ2((/mach1);

Fig.11 Profil aerodinamic aripa

Curgere trans-sonica

II.4.6 CONVECȚIE

Fie u(x) un câmp de viteze iar ø concentrația unui produs chimic. produsul chimic poate fie să fie prezent în momentul zero sau să fie produs în interiorul fluidului de către o sursă f sau să fie generat la frontieră de către o sursă øΓ sau o curgere g.

ecuația pentru evoluția temporală și spațială a ø este

în Ωx]0, T [

unde v este coeficientul de difuzie. aici fluidul este incompresibil, altfel u· ar fi înlocuit de .

următoarele condiții de frontieră vor defini în mod unic ø

ø(x, 0 ) = ø0 ( x )

Γ1 U Γ2 = Γ Γ1 ∩ Γ2 = 0

Sunt posibile câteva scheme de discretizare a timpului; schema Crank – Nicolson prezintă precizie și stabilitate:

Am presupus că u și f sunt independente de timp. introducând variabila auxiliară ω = (øm+1 + øm ) / 2 putem scrie schema ca

;

Aceasta s-ar programa ca și:

iter ( … )

begin …

solve(w) id(w)*2/dt+dx(w)*u1+dy(w)*u2-laplace(w)*nu=f+phi*2/dt;

phi = 2*w-phi;

end

Când v este prea mic metoda produce oscilații de aceea este problematic.

Înainte de toate, există cazuri unde v=0 este fizic. de exemplu, uleiul nu se amestecă cu apa și astfel, peste domenii mici, termenul de difuziune nu există când se studiază transportul și dispersia uleiului de către un curent de apă. Ecuația este atunci cunoscută drept convecție, sau transport și este bine poziționată cu următoarele condiții de frontieră:

în Ωx]0, T [

ø(x, 0 ) = ø0 ( x )

Γ1 U Γ2 = Γ – Γ1 ∩ Γ2 = 0

unde Γ – este parte a frontierei unde curgerea intră Ω

Γ –= { xЄΓ: u(x).n(x) < 0 }

Acum, aceasta sugerează în mod puternic, că trebuie luată informația din amonte în schemă ( așa numită upwinding ). o modalitate este de a aminti că

la t = tm+1

pentru programare există un operator special care transformă orice funcție ø în ø ( x – u (x) dt ); este chemată de către convect (phi, u1, u2, dt ). Deci, o altă schemă pentru problemă este

iter ( … )

begin…

solve(w) id(w) * 2/dt -laplace(w) * nu

= f + phi/dt + convect(phi,ul,u2,dt)*/dt;

phi = 2 * w – phi;

end

Pentru ilustrare luați în considerare o porțiune de conductă și luați în considerare o curgere potențială. Conducta este duza studiată mai sus cu excepția că este separată în două ramuri. o concentrare inițială de produs chimic este zero peste tot cu excepția lângă punctul 0.8a, 0. Curgerea poartă produsul chimic iar separarea în două ramuri poate fi studiată:

a:= 6; b:= 1; c:=0.5; pi := 4*atan{l):

changewait;

border(l,0,l,8)begin x:=-a; y:=l+b–2*(l + b)*t end;

border(2,0,1,26) begin x: = -a+2*a*t;

y:=-1-b*(x/a)*(x/a)*(3-2*abs(x)/a)

end;

border(3,0,1,8)begin x:=a;y:=-l-b+(1+ b)*t end;

border(4,0,1,20) begin x:= a – a*t; y:=0 end;

border(4,0,pi,8)begin x:=-c*sin(t)/2;y:=c/2-c*cos(t)/2
end;

border(4,0,1,30) begin X:= a*t; y:=c end;

border(3.0,1,8)begin x:=a;y:=c +(1+ b-c)*t end;

border(5,0,1,35) begin x:= a-2*a*t;

y:=l+b*(x/a)*(x/a)*(3-2*abs(x)/a) end;

buildmesh(800);

onbdy(1) psi = y;

onbdy(2) psi =-1 -b;

onbdy(3) dnu(psi) =0;

onbdy(4)psi=0; /*controlează cum se separă curgerea*/

onbdy(5) psi =l+b;

solve(psi) -laplace{psi) =0;

plot(psi);

ul = dy(psi); u2 = -dx(psi); dt =0.3; nu=0.01;

onbdy(1)w=0;

onbdy (2)w=0;

onbdy(3)dnu(w)=0;

onbdy(4)w=0;

onbdy(5)w=0;

solve(w,l)id(w)*2/dt-laplace(w)*nu=exp(-10*((x+a*0.8)ˆ2+(yˆ2)));

phi =w;

plot(phi);

iter(15)

begin

f = convect(phi,ul,u2,dt);

plot (f) ;

solve(w,-l) id(w) * 2/dt -laplace(w) * nu

= phi/dt +f/dt;

phi = 2 * w – phi;

plot(phi);

end

Trebuie specificat că nu am repetat condițiile de frontieră în interiorul buclei pentru ω. din moment ce matricea este deja factorizată, vor fi același tip de condiții de frontieră cu zero pe dreapta, ceea ce ne dorim. Acest lucru se întâmplă pentru că toate funcțiile sunt inițializate la zero

printr-un interpret

Fig.12 Fenomenul de curgere-separarea

II.4.7 ECUAȚII NAVIER – STOKES

7.1. Generalități

Curgerile vâscoase incompresibile satisfac :

în Ω]0,T[

u׀t=0=u0, u׀Γ = uΓ

Un algoritm posibil este

,

unde uoX(x) = u (x – u(x)δt). Unele modificări trebuie să fie incluse când frontiera aproximează infinitul. Atunci, dacă este o frontieră de aflux, u, ∂np sunt date și dacă este o frontieră de eflux ∂np, p pot fi date, dacă e posibil.

Problema de pas înapoi

Metoda a fost testată pentru curgerea peste un pas. Acesta funcționează când triangulația este suficient de bună pentru un anumit număr Reynolds dat.

Observație :

Din moment ce interpretul netezește derivatele în funcții P1 când dx(u) este invocat, trebuie avut grijă pentru a reduce difuzia numerică. de exemplu, pentru această problemă, calculul lui pare a fi prea răspândit / difuzabil pentru a calcula funcția de curent, ψ corect prin rezolvarea

-Δψ=dy(u)-dx(v)

changewait;

border(1,0,1,6) begin x:=0; y:=l-t end;

border(2,0,1,15) begin x:=2*t; y:=0 end;

border(2,0,1,10) begin x:=2; y:=-t end;

border[2,0,1,20) begin x:=2+3*t; y:=-l end;

border(2,0,1,35) begin x:=5+15*t; y:=-l end;

border(3,0,1,10) begin x:=20; y:=-l+2*t end;

border(4,0,1,35) begin x:=5+15*(1-t); y:=l end;

border(4,0,1,40) begin x:=5*{l-t); y:=l end;
buildmesh(900);

nu = 0.002; dt := 0.4;

/* presiune inițială */

onbdy(1)dnu(p) =-2*nu;

onbdy(2) dnu(p) = 0;

onbdy(3) p=0;

onbdy{4) dnu(p) = 0;

solve(p,l) – laplace(p)= 0;

/* viteză inițială orizontală */

onbdy(1) u = y*(1-y);

onbdy(2) u = 0;

onbdy(3) dnu(u) = 0;

onbdy(4) u = 0;

solve(u,2)id(u)/dt-laplace(u)*nu=-min(y*y-y,0)/dt;

/* viteză inițială verticală */

onbdy(1)v = 0;

onbdy(2) v = 0;

onbdy(3) v = 0;

onbdy(4) v = 0:

solve(v,3) id(v)/dt-laplace(v)*nu =0;

un = u; vn = v;

iter(80)

begin

f=convect(un,u,v,dt); g=convect(vn,u,v,dt);

/* Viteză orizontală */

onbdy(1) u = y*(l-y);

onbdy(2) u = 0;

onbdy(3)dnu(u)=0;

onbdy(4) u = 0;

solve(u,-2) id(u)/dt-laplace(u)*nu = f/dt –dxt(p);

plot(u);

/* Viteză verticală */

onbdy(1) v = 0;

onbdy(2) v = 0;

onbdyd) v = 0;

onbdy{4) v = 0;

solve(v,-3) id(v)/dt-laplace(v)*nu=g/dt–dy(p);

/* Presiunea */

onbdy(1)dnu(p) a-2*nu;

onbdy(2) dnu(p) = 0;

onbdy(3) p=0;

onbdy(4) dnu(p) = 0;

solve(p,-l) -laplace(p)= -(dx(f) + dy(g))/dt;

un = u;

vn = v;

end ;

save('u.dta',u);save('v.dta',v);save('p.dta',p)

plot(u);

Fig.13 Liniile de nivel ale vitezei orizontale pentru o curgere Navier – Stokes prin intermediul unui backward step Re=100, ( bazat pe viteza la centrul frontierei de admisie și înălțimea pasului ).

7.2. Straturile frontierei

Din moment ce curgerea dezvoltă gradienți puternici lângă frontieră, mesh-ul trebuie să fie îmbunătățit în direcția normală lângă frontieră. acest lucru este efectuat prin adăugarea de frontiere artificiale lângă cele fizice. Aici am adăugat una ușor deasupra frontierei de jos. este important a se atribui un număr 0 acestei frontiere altfel devine una fizică.

/* Frontiera fizică */

border(1,0,1,6) begin x:=0; y:=l-t end;

border(2,0,1,15) begin x:=2*t; y:=0 end;

border(2,0,1,10) begin x: = 2; y:= -t end;

border(2,0,1,20) begin x: = 2 + 3*t; y:=-1 end;

border(2,0,1,35) begin x:=5+l5*t; y:=-l end;

border[3,0,1,10) begin x:=20; y:=-l+2*t end;

border(4,0,1,35) begin x:=5+15*(1-t); y:=l end;

border(4,0,1,40) begin x:=5*(1-t) ;y:=1 end;

/* O frontieră artificială deasupra de bdy 2 la distanța d */

d := 0.03;

border(0,d,l+d/2,15)begin x:=2*t; y:=d end;

border (0,0,1, 10)begin x:=2+d; y:=d-t end;

border(0,0,1-d,20)begin x:=2+d*3*t;y:=-l+d end;

border(0,0,l-d/5,35)begin x:=5-2*d+15*t,y:=-l+d end;

/* Alta la distanța de 4*d */

d := 4*d;

border(0,d,l+d/2,15)begin x:=2*t;y:=d end;

border (0,0,1,10)begin x:=2+d; y:=d-t end;

border(0,0,1-d,20)begin x:=2+d+3*t;y:=-l+d end;

border(0,0,1-d,35)begin x:=5-2*d+15*t;y:=-l+d end;

baildmeshd(1500)

Fig.14 Adaugarea unei frontiere artificiale

CAPITOLUL III Modelarea Fenomenelor

de camp electromagnetic

Cuprins

III.1 Introducere

III.2 Metode integrale și descompunere modală

III.3 Metoda elementului finit

III.4 Modelul de cavitate cu corp de proba

III.5 Determinarea profilelor de camp electromagnetic

III.1 Introducere

În probleme de câmpuri electromagnetice, cuplate cu cele termice, descrierea procesului se face cu ajutorul ecuațiilor diferențiale, a căror rezolvare este dificilă din cauza complexității structurilor de studiat. Datorită dezvoltării metodelor de calcul numeric, s-au pus la punct mai multe moduri de abordare a determinării caracteristicilor și a rezolvării ecuațiilor diferențiale.

Trebuie recunoscut interesul pe care îl prezintă metodele variaționale și metodele integrate, care sunt frecvent folosite în modelarea circuitelor cu microunde, permițând abordarea problemelor de limite și discontinuități.

O altă metodă, cea a elementului finit, oferă rezultate concludente în rezolvarea ecuațiilor de câmp și a ecuațiilor de căldură. Aceasta permite caracterizarea sistemelor complexe, cum ar fi cavități rezonante, cu sarcini plasate aleator, oricare ar fi forma acestora.

In figura A, se prezinta legatura existenta intre metodele de calcul numeric.

Fig A Relații între metodele numerice existente.

III.2 Metode integrale și descompunere modală

Aceste metode sunt utilizate pentru a modela tranziția ghid- cavitate și unele discontinuități metalice.

a. Metoda curenților echivalenți

Această metodă are ca principiu, înlocuirea discontinuităților metalice cu densitatea de curent echivalentă, care se comportă la fel ca o antenă în interiorul structurii.

Fig. III.1. Difracția câmpului, în prezența unui obstacol metalic

Pornind de la ecuația de propagare a câmpului electric sub forma :

În formula de mai sus, nu mai este necesar să explicităm vectorii asupra cărora se aplică operatorii, pentru că toate câmpurile sunt soluția aceleași ecuații [5]. Ecuația (III-l), se poate rezolva cu ajutorul funcțiilor lui Green diadice, care satisfac următoarea ecuație :

unde (M- Mo)-funcția lui Dirac, care este nulă în afară de M=Mo

I- diadic vector unitate

Funcția lui Green diadică este obținută prin rezolvarea ecuației (III-2), prin descompunerea sa în funcții proprii ale structurii respective:

Câmpul electric radiat de curentul de pe discontinuitate este:

Cîmpul total în structură se obține prin suprapunerea câmpului de excitație Einc

cu câmpul difractat de discontinuitate

Prin înlocuirea (III-3), în (III-4), se poate obține câmpul total, descompunându-l într-o serie de moduri ale ghidului:

Pornind de la expresiile (III-6) și (III-7), într-un ghid monomod se poate determina valoarea impedanței acestuia.

b. Metoda integralei de frontieră

Vom folosi această metodă pentru a modela joncțiunea ghid–cavitate. Cu ajutorul acestei metode, problemele de volum se pot reduce la probleme de frontieră. Considerăm o cavitate de volum v, ca în Fig.III.2 limitată de o suprafață conductoare V și excitată prin deschiderile Si

Fig. III.2. Cavitate, excitată cu mai multe surse

Scriind ecuațiile lui Maxwell pentru volumul v avem:

Ținând cont de condițiile la limită din punct de vedere electric și magnetic:

Ecuația de propagare ce rezultă din ecuațiile lui Maxwell, devine :

Funcția lui Green diadică, verifică ecuația:

Înmulțind apoi, relația (III-10), la dreapta scalar, cu G și relația (III-11) la stânga scalar cu E și scăzând membru cu membru obținem:

și ținând cont de identitatea lui Green diadică:

Integrând această ultimă ecuație și aplicând apoi formula lui Ostrogradsky [126], obținem:

unde: S = S1+ S2 + S3

iar

Funcția lui Green diadică, trebuie să respecte condițiile la limită impuse pe suprafața v .Dacă aceste condiții sunt de tip electric, atunci funcția lui Green este notată , dacă sunt de tip magnetic notația este . Ecuațiile (III-14), (III-15), devin :

Cu ajutorul acestor expresii, cunoscând câmpul pe deschiderile Si, putem determina câmpul electromagnetic în orice punct al volumului v [67]. Funcția lui Green diadică, este soluția ecuației (III-1o), ținând cont de condițiile la limită, impuse de conturul volumului v. Tocmai aceste condiții la limită permit determinarea funcțiilor proprii ale cavității, numite moduri. Acestea formează o bază de descompunere pentru funcția lui Green, asemănătoare cu cea dată de (III-3):

Înlocuind (III-18) cu (III-11), amplitudinile se obțin prin aplicarea condițiilor de ortogonalitate proprii structurii, rezultă:

Caractreristica principală a structurii prezentată în Fig.III.2. este o matrice de difracție, care trebuie determinată pornind de la ecuațiile (III-19), (III-2o). Pentru aceasta este necesară scrierea ecuațiilor de câmp ca o combinație liniară de moduri proprii, relativ la fiecare ghid de excitație:

Considerând, pentru exemplificare, condițiile electice, asupra unei cavități, vom ține cont că pe suprafața Si, câmpul electric tangențial se anulează și deci integrala (III-16), este zero. Identificarea se face cu ajutorul integralei (III-17), utilizând funcția lui Green Gh, din (III-2o),

unde l = 1,2,3

La nivelul fiecărei joncțiuni ghid – cavitate avem egalitate între (III-23) și (III-22). Proiectând fiecare egalitate asupra modurilor ghidurilor de acces respective, putem deduce matricea de difracție a multipolului [77].

B = ( S ) A (III-24)

unde A matricea coloană a modurilor incidente spre cavitate Aim

B matricea coloană a modurilor reflectate de către cavitate Bim

( S ) matricea de difracție

Fig. III.3. Joncțiunea ghid-cavitate “ i “

Desigur că, diferitele serii care intervin în calcul se elimină în funcție de convergența elementelor matricei.

III.3 Metoda elementului finit

Această metodă, ca de altfel și alte metode numerice , necesită o formulare variațională, adică trebuie căutată minimizarea unei funcții ce rezultă din ecuațiile diferențiale care descriu sistemul. Iată câteva tipuri de formulări:

a- Formularea variațională

Având de rezolvat o ecuație de tipul

unde L – este un operator cunoscut

g -expresie cunoscută

f -funcția necunoscută

Funcția (III-25), reprezintă ecuația de undă (III-1). Soluția problemei este o funcție f, care face să fie staționară, în raport cu mici variații , o expresie de formă:

care se poate scrie sub formă integrală :

unde F =f*Lf-(f*g+fg*), ( )* – complex conjugat.

Cum integrandul lui F ( f ) este staționar, înseamnă că pentru orice variație f, avem

De îndată ce putem găsi un principiu variațional, vom putea obține imediat o aproximare sub formă integrală.

Considerând o funcție n , de forma:

Datorită faptului că I ( f ) este staționar avem în consecință :

Deoarece relația (III-29), trebuie verificată pentru orice variație an, vom obține un sistem de ecuații, care poate duce la rezolvarea necunoscutelor an și prin urmare a funcției aproximative fapr

b- Formularea proiectivă

O altă metodă, care nu apelează la formularea variațională, este metoda proiectivă, numită și metoda restului ponderat. Pornind de la realația (III-25), care se poate scrie și sub forma :

Considerând fapr, soluția aproximativă a acestei ecuații, restul se poate scrie :

Metoda se bazează pe o treoremă specifică spațiilor Hilbert, care stabilește că într-un asemenea spațiu, numai vectorul nul este ortogonal cu ceilalți vectori ai spațiului. Bazându-ne pe acest principiu o soluție de forma (III-28) este perfect posibilă, fiind suficient să facem proiecția restului pe o bază de funcții ponderate n

Contrar metodei varianționale, metoda restului ponderat , nu conduce întotdeauna la matrici simetrice. Rezolvarea numerică, pornind la aceste două metode este dificilă din cauza mărimii matricii care modelează ecuația diferențială.

Aceste modele de ecuații se mai pot rezolva și cu metoda lui Galerkin, care este un caz particular al metodei proiective.

c-Formularea celor mai mici pătrate

Pentru rezolvarea unor probleme fizice, având condiții suplimentare , metoda celor mai mici pătrate, reprezintă o alternativă viabilă, pentru obținerea unor formulări integrale, prin care să se poată construi o aproximație.

Această metodă are ca principiu, minimizarea sumei pătratelor resturilor ecuației diferențiale. Pornind de la ecuație diferențială (II-30), se definește mărimea:

Căutând să obținem o soluție de forma (III-28), vom găsi o rezolvare corectă, minimizând expresia de forma:

Matricea care decurge din această formulare este simetrică și pozitiv definită.

d-Rezolvarea cu ajutorul elementului finit

Rezolvarea unui sistem diferențial de forma (III-25), are ca primă etapă divizarea domeniului , în subdomenii e , numite elemente. Metoda se bazează pe formulele integrale amintite, utilizând ca funcții de bază, funcțiile de frontieră. Datorită acestei noi funcții , se face rezolvarea unei ecuații matriciale ale cărei matrici prezintă elemente nule, numită matrice încrucișată[30].

Funcțiile de aproximare sunt polinoame de interpolare, iar valoarea funcției necunoscute pe modurile unui element o numim valoarea modală.

Reluând ecuația diferențială (III-25), care este însoțită de anumite condiții de limită ( CL ):

În problemele de termica microundelor L reprezintă ecuația de undă a câmpului sau ecuația căldurii, iar CL reprezintă condițiile la limită de tip Dirichlet sau Newman .

Vom aplica formularea integrală ( III-31 ), în metoda elementului finit, pentru fiecare subdomeniu e :

N matricea linie a funcțiilor de interpolare

{fn }si {n } matrice coloană a variabilelor nodale

Deoarece matricea N , este zero în exteriorul elementului e și fn nu face să intervină decât variabilele nodale ale aceluiași element, fiecare termen Ie se calculează pornind de la variabilele legate de acest element.

Datorită naturii repetitive a acestor operații, proprietatea descrisă, a contribuit la ușurarea calculului numeric și la promovarea cu succes a metodei elementului finit.

Rezultatul global se obține, datorită continuității funcției f, prin :

unde Nt numărul total de subdomenii.

La fel ca integrala I, dată de (III-27), gruparea tutor elementelor e, a domeniului , permit obținerea unei ecuații matriciale, ale cărei necunoscute sunt variabile nodale

unde L matrice globală

{f} vector global al tuturor variabilelor nodale

{g} vector global al tuturor condițiilor impuse

În consecință pentru găsirea soluției finale, trebuie să rezolvăm ecuația matricială ( III- 35 ).n ceea ce privește convergența soluției aproximative, spre soluția problemei fizice, există câteva mijloace pentru diminuarea erorilor:

Creșterea numărului de subdomenii și micșorarea dimensiunilor lui e Soluția convergentă se obține atunci când mărimea lui e tinde la zero.

Mărirea gradului de libertate al polinoamelor care reprezintă funcțiile de interpolare. Eroarea tinde la zero când gradul de libertate tinde la infinit.

III.4 Modelul de cavitate cu corp de proba

Vom considera un cuptor de microunde în care se află o sarcină de formă oarecare și vom încerca să aflăm distribuția puterii disipate în sarcină, prin două metode: Prima, o metodă hibridă, care folosește descompuneriea modală combinată cu metoda elementului finit, iar a doua, metoda elementului finit singură.

metoda descompunerii modale

metoda descompunerii modale

metoda elementului finit

Fig. III.4. Cavitate cu sarcină disipativă

Metoda descompunerii modale este utilizată în volumul V2 , care nu conține sarcina, iar metoda elementului finit este folosită în volumul V1 care conține sarcina, deoarece această metodă se adaptează foarte bine la orice formă geometrică. La separația dintre volumele V1 și V2 se crează o suprafață fictivă, care face racordul dintre câmpuri și permite determinarea unei matrici de difracție.

În Fig. III.4., avem o cavitate rezonantă, în care se află un mediu disipativ, caracterizat prin ’ , ” .Vom face o divizare fictivă a cavitații în volumele V1 și V2, separate printr-o suprafață fictivă S1. Metoda hibridă pe care o vom folosi, va duce la înlocuirea volumului V1, cu o matrice de difracție echivalentă .

Deoarece metoda va fi utilizată într-o aproximare 2D, atât structura pe care o vom studia, cât și excitația sunt independente de axa y. Câmpul electric are o singură componentă, anume Ey:

câmpul electric este soluția ecuației

iar Ey = 0 pe pereții metalici.

Vom afla mai întîi câmpul pe suprafața S1, pentru ca apoi cu metoda elementului finit, să aflăm câmpurile în volumul V1.

Pentru aceasta în volumul V2, descompunem câmpul într-o serie de moduri TEno . Astfel pentru un mod incident de ordin m, câmpurile tangențiale în vecinătatea S1 se pot scrie:

Admitanța Y1n este reală dacă modul se propagă și imaginară dacă este la tăiere. 1n este constanta de propagare a modului n la interfața S1, condițiile de excitare vor fi impuse sub forma:

unde N1-numărul de moduri existente la nivelul S1.

Cu ajutorul elementului finit, se calculează E1,n în interiorul volumului V1, iar apoi se deduce câmpul magnetic , printr-o derivare numerică.

Ne interesează numai componeneta tangențială a câmpului magnetic, care are următoarea formă:

Pornind de la câmpul electric pe S1, dat de (III-38), prin metoda elementului finit vom ajunge la câmpul magnetic:

Ținând cont de continuitatea câmpului magnetic de o parte și cealaltă a lui S1, egalând expresiile (II-39) și(III-45), identificând termenii, obținem:

Rezolvarea sistemului liniar, duce la obținerea liniei m, a matricii de difracție {S1mn} . Trebuie remarcat că dacă suprafața de separție S1, este depărtată de sarcină, atunci este suficient să luăm în considerare modurile care se propagă în cavitate, altfel trebuie ținut cont și de primele moduri de tăiere.

Ordinul lui{S1mn}, depinde de numărul de moduri existente la nivelul suprafeței S1. [30]. Putem în continuare, să înlocuim volumul V1, cu matricea de difracție echivalentă:

unde {A1n} – matrice coloană a amplitudinilor de câmp care ies din V1 la nivelul S1.

{B1n} – matrice coloană a amplitudinilor de câmp care intră în V1, la nivelul S1

Fig. III.5. Cavitatea cu matricea de difracție

Volumul V1 este înlocuit de matricea de difracție {S1mn}, la nivelul suprafeței S1. Câmpul în cavitate este funcție de excitația din ghidul de undă, el descompunându-se într-o serie de moduri. Pentru studiul structurii din Fig.II.5., folosim metoda integralei de frontieră [8], care permite calculul câmpului în orice punct, în funcție de câmpul la cele două intrări ale cavității.

Reluând expresiile (III-16) și(III-17), le putem aranja sub forma:

S0=S1+S2 suprafețele unde cîmpul este diferit de zero.

La nivelul celor două accese, cîmpurile electrice și magnetice se pot scrie, după cum urmează [5],[66]

Alegerea condițiilor la limită asupra funcției lui Green, permit eliminarea uneia dintre integralele ( III-48 ) sau( III-49 ). Impunând condiția ca pe S0 câmpul electric să fie nul, rezolvarea se bazează pe sistemul dat de ( III-49 ) și

( III.51 ), cu funcția lui Green diadică, rezultată din ( III-20 ) :

Funcția lui Green, verifică următoarea condiție:

și deci integrandul din (III-49), se poate scrie:

Se poate observa că, în fiecare acces, ținem cont doar de câmpurile tangențiale, componentele normale ale câmpului nu intervin.

Pentru găsirea soluțiilor, va trebui să rezolvăm sistemul de ecuații obținut folosindu-ne de continuitatea câmpului magnetic și de deschiderile S1 și S2:

unde l= 1,2-desemnează suprafețele S1 și S2

i=1,2-desemnează accesele

Suntem în fața a două ecuații (i=1,2), fiecare conținând două integrale (l=1,2), care trebuiesc calculate. Pentru simplificarea calculului numeric, vom încerca să suprimăm diferite serii.

Pentru i=1, vom obține:

Vom înlocui acum, H și Gh cu expresiile lor și ne vom folosi de câțiva operatori pentru a transforma egalitatea (III-58) în (III-59) :

Pentru a determina relațiile între amplitudinea Alj și Blj, este suficient să proiectăm relația (III-59), pe baza de funcții proprii ale accesului excitației Vom obține N1 ecuații, ale căror necunoscute sunt amplitudinea câmpurilor. Utilizând apoi proprietatea de ortogonalitate, vom avea, pentru j=1,N1

O relație similară cu (III-60), se poate obține, de la sistemul (III-56) și (III-57), pentru i = 2, ceea ce ne permite să adăugăm N2 ecuații, la cele N1 deduse din relația (III-60).

Sistemul acesta poate fi pus sub formă matricială, ca o expresie între amplitudinile câmpurilor care ies, în funcție de cele care intră:

Pentru a exprima o soluție generală a problemei și pentru a face legătura între volumele V1 și V2 ale cavității, trebuiesc combinate cele două relații matriciale (III-47) și (III-61), care pot fi puse sub forma:

Din aceste două relații matriciale, se poate deduce relația care modelează întreaga cavitate:

Deoarece cavitatea, este în general excitată printr-un singur mod fundamental al ghidului de acces, matricea coloană {A2m}, are doar primul element diferit de zero.

Pentru ca acum să calculăm câmpul în volumul V1, va trebui să procedăm invers. Cunoscând excitația, în cazul de față A12, sistemul (III-62), permite determinarea amplitudinilor de câmp, la nivelul suprafeței S1, care sunt:

Înlocuind apoi, aceste amplitudini cu valorile lor în (III-50), vom obține valorile de câmp la suprafața S1. În continuare, cu metoda elementului finit, putem obține valorile de câmp în orice punct al volumului V1 [44].

Compararea metodei de lucru

Pentru a putea face o comparație între două metode de abordare, una hibridă (metoda elmentului finit, combinată cu metoda elementului de frontieră), iar cealaltă, metoda elementului finit, vom considera următoarea aplicație:

Este vorba despre o cavitate, în care se află un mediu disipativ, așa cum se vede în Fig.III.6.

Fig. III.6. Cavitate reală cu dielectric

În Fig.III.7., am reprezentat , distribuția puterii disipate în sarcină, la o adâncime de 9 mm, de la suprafața (), rezultate obținute prin metoda hibridă și prin metoda elementului finit [54].

Fig. III.7. Distribuția puterii disipate în sarcină la 9mm de suprafață ()

Se poate constata că diferențele între rezultatele obținute prin cele două metode, se situează sub 5,3%. Pentru metoda hibridă, la descompunerea modală s-au folosit 11 moduri la joncțiunea cavitate – ghid, având ca scop convergența rezultatului, iar pentru metoda elementului finit, s-au folosit 3900 elemente.

III.5 Determinarea profilelor de camp electromagnetic

Aplicarea programului Mac-GFEM

Metoda elementului finit, este bazată pe subdivizarea unui domeniu , într-un număr finit de elemente e . În cazul softului Mac-GFEM, am utilizat un generator de rețea 2D, care folosește o interfață grafică Macintosh, pentru a crea rețele cu ajutorul triunghiurilor. Rețeaua este înregistrată în text neformatat, ceea ce permite reducerea ocupării memoriei. [146 ]

Etapele care trebuiesc urmate sunt:

-Definirea frontierelor structurii

-Divizarea în triunghiuri și generarea rețelei

-Punerea condițiilor la limită

-Definirea ecuațiilor de rezolvat

-Vizualuizarea rezultatelor

Pentru a putea simplifica definirea analitică a elementelor de formă complexă, vom introduce noțiunea de element de referință. Acesta se definește ca un element de formă simplă r, situat într-un spațiu de referință și care poate fi transformat în oricare element real e, printr-o transformare geometrică Ze

Astfel în cazul unui triunghi:

Transformările e, trebuie să genereze elemente reale, care să satisfacă criteriul următor :

De aceea , fiecare transformare trebuie să aibă următoarele proprietăți :

să fie bijectivă în orice punct situat în interiorul elementului de

referință sau pe frontiera sa

nodurile geometrice ale elementului de referință, corespund celor ale elementului real

fiecare porțiune de frontieră a lui r, definită de nodurile geometrice respective, corespund porțiunii de frontieră a lui e, definit de nodurile corespunzătoare.

Transpunerea e este bijectivă, dacă matricea jacobiană J nu este singulară, adică det (J) este diferit de zero,

Determinantul acestei matrici nu se poate anula, decît atunci când cele trei noduri geometrice sunt coliniare.

Considerând acum funcția f(x), vom putea să o aproximăm convenabil în felul următor:

Vom diviza n intervale, (b-a)/n= și luănd pe fiecare interval o funcție constantă, se obține noua funcție aproximativă fa ;dacă 0 , atunci fa(x)f(x)

Dacă în fiecare interval se lucrează cu funcții binome ax+b, atunci vom obține :

Dacă dorim ca funcția aproximativă să conveargă mai repede, trebuie să folosim polinoame cu grad mare.

Lucrând în două dimensiuni, vom proceda la divizarea planului în mici triunghiuri ca în reprezentarea următoare :

Dacă se cunoaște funcția pe noduri,atunci se poate interpola pe triunghiuri [38].

Modelarea cavității cu sarcină folosind metoda elementului finit în 2D

Din punct de vedere al câmpului electromagnetic, considerând ecuația de propagare :

folosind formularea integrală pe domeniul , divizat în subdomenii i, vom obține:

Pentru un triunghi i

Vom ajunge astfel , la o ecuație de forma:

AX=B

unde X este o matrice coloană de câmpuri pe noduri. Vom putea determina astfel valorile lui E în noduri și deci profilul electric în spațiul dorit.

Pentru domeniul în care vrem să obținem profilul termic, ecuația căldurii este [110]:

Pe frontierea d, condițiile la limită sunt:

Presupunând Tn cunoscut în (x,y), prin interații pornind de la condițiile inițiale, ne rămâne de rezolvat ecuația:

deoarece într-un mediu eterogen, a este funcție de coordonate: a(x,y).

Conform metodologiei elementului finit, se aproximează funcția T, în fiecare element definit de rețea:

unde i este funcția de bază

Acestă funcție i , este de ordinul întîi, adică se presupun în fiecare element, variații liniare ale lui T. Deasemenea funcția de bază este nulă în celelalte elemente. În cele N elemente ale lui , funcția T este :

Metoda de rezolvare, constă în căutarea coeficienților bi, proiectând operatorul pe funcțiile de bază, conform metodei lui Galerkin.

Algoritmul elementului finit, aduce în final ecuația căldurii la o expresie matricială, care permite determinarea temperaturilor aparținând fiecărui nod din rețea, ținând cont de condițiile de limită și de surse. Cunoscând apoi temperaturile pe noduri, se pot calcula temperaturile în fiecare element, presupunând că variațiile sunt liniare și de aici vom putea determina profilul termic în întreg domeniul.

Matricea obținută, de dimensiuni NxN, este o matrice de bandă, simetrică, ceea ce are avantajul limitării timpului de calcul [110].

Determinarea profilelor câmpului electromagnetic

Așa cum se va vedea în Cap.V., în urma calculelor de dimensionare a cavității rezonante, ținînd cont de geometria creuzetului ceramic (transformator de căldură) și de corpul de probă ( sarcină absorbantă), s-a ajuns la dimensiunile geometrice din Fig.III.8. și la amplasarea din Fig.III.9.

Programul utilizat fiind MacGFEM în 2D, vom parcurge etapele descrise în Cap.III.5.1.

Odată ce frontierele structurii au fost definite, vom trece la divizarea în triunghiuri și generarea rețelei.

Fiecare contur, denumit în continuare regiune, a fost divizat în segmente în funcție de lungimea acestora și de densitatea rețelei pe care am dorit să o obținem [11],[42]

Programul Mac-GFEM, poate lucra cu maximum 4000 de noduri și 10 regiuni, ceea ce a adus la necesitatea unei realizări foarte atente a rețelei de triunghiuri.

Faptul că generatorul de rețea alege automat un număr mai mare de triunghiuri pentru segmentele scurte ale oricărui contur, a impus necesitatea unei segmentări judicioase a cavității, a transformatorului de căldură și a sarcinii, astfel încât să avem o repartiție uniformă a rețelei în întreaga cavitate și de asemenea o concentrare de triunghiuri în zonele care prezintă interes pentru vizualizare.

În Fig.III-10., se poate vedea că în zona transformatorului de căldură și a sarcinii există un maxim de triunghiuri, astfel încât informațiile privind câmpul electromagnetic și temperaturile vor fi concentrate în această zonă.

Trecând acum la condițiile la limită, am considerat că izolația termică dintre transformatorul de căldură și pereții cutiei metalice este suficient de bună, astfel încât să nu avem pierderi de căldură, temperatura va fi deci constantă, pe suprafața exterioară a transformatorului de căldură.

Fig. III.8. Cavitate rezonantă, cu transformator de căldură și sarcină absorbantă. Notații geometrice.

Fig. III.9. Cavitate cu transformator de căldură și corp de probă.

Planul de simetrie

Se poate observa în Fig.III.8., că transformatorul de căldură are o latură deschisă pentru a putea introduce corpul de probă. Ansamblul va fi înconjurat de un capac virtual care să ne dea posibilitate să punem condițiile la limită, respectiv temperatura Tp = const. Acest capac virtual va avea proprietățile dielectrice ale aerului, astfel, nu va concentra câmp și nu va acumula căldură, putând fi considerat un corp neutru electromagnetic și termic.

Prin ghidurile de undă, dispuse central, se vor propaga undele electromagnetice, joncțiunea cu cavitatea făcându-se în plane diferite, având o poziționare care să evite fenomenul de cuplaj electromagnetic între surse.

Tinând cont de condițiile de frontieră impuse, câmpul electromagnetic se va concentra în transformatorul de căldură și sarcină, iar transferul de căldură va avea loc în mediul umplut cu aer dintre cele două.

Ecuațiile care trebuiesc rezolvate sunt:

– pentru ghidul 2

unde am considerat că u este echivalent cu H , u = 1 pe frontiera (2), pe toate celelalte frontiere ale incintei, acesta fiind nul.

r= permitivitatea relativă a mediului

c = viteza luminii

f = frecvența- 2,45 l09 Hz

– pentru ghidul 1

am considerat că v este echivalent cu E

pe frontiera (4)

pe toate celelalte frontiere ale incintei acesta fiind nul.

Pentru a putea face superpoziția celor două câmpuri, vom face transformarea între H și E astfel:

Deoarece ghidurile notate Magl și Mag 2 sunt identice, vom căuta relația dintre E0 și H0, astfel :

Deasemenea considerând neglijabile pierderile în vid, și ținând cont că generatoarele de microunde au puterea maximă de l kW, vom putea obține câmpul electric într-o secțiune prin cavitate astfel:

Valorile lui r1 si r2., au fost preluate în urma determinărilor experimentale pe eșantioane de ceramici, conform Cap.IV.

Explicația pentru care am optat să lucrăm pentru ghidul l, cu câmpul electric E iar pentru ghidul 2, cu câmpul magnetic H, ține de posibilitățile de lucru ale softului utilizat.

Astfel Mac-GFEM, poate lucra numai cu funcții continue, ori așa cum se vede din Fig.III.9., câmpul Ey din ghidul 2, oscilează într-un același plan cu sarcina, în secțiunea transversală prin structură, ceea ce duce la apariția unei discontinuități. Pentru a evita această situație, vom lucra cu câmpul magnetic H în ghidul 2, a cărui orientare este perpendiculară pe planul sarcinii, putând să-l tratăm ca pe o funcție continuă.

În situația în care cele două magnetroane vor funcționa simultan câmpul total în fiecare punct din incintă, va fi dat de superpoziția lui E și H , ținând cont că sunt funcții complexe, obținem mărimea Etot.

Vom ține cont mai departe că ” r = Im(r )

unde r- reprezintă o funcție ce determină valoarea permitivității

în orice regiune a structurii ,astfel :

r=1 -pentru regiunile o,3 și 4 umplute cu aer

r =r1 -pentru regiunea cu sarcina

r =r2 -pentru regiunea cu transformatorul de căldură

Ajungem , folosind relația (II-21), în aproximația că pierderile în dielectrici sunt date numai de partea imaginară a lui r , la expresia puterii active de microunde, disipată într-un material, într-un punct M, la momentul t,

Vom putea obține acum profilele câmpului electromagnetic, pe care le vom prezenta comparativ, cavitate cu sarcină, cu sau fără transformator, cu unul sau două generatoare de microunde, în Fig.III.10, Fig.III.11, Fig.III.12, Fig.III.13, Fig.III.14, și Fig.III.15. prezentate în paginile ce urmează :

Fig.III.10 Mesh-reteaua de triangulatie in incinta

Fig.III.11 Distributia campului electromagnetic

Functionare magnetron sus

Fig.III.12 Distributia puterii disipate

Functionare magnetron sus

Fig.III.13 Distributia campului electromagnetic

Functionare magnetron dreapta

Fig.III.14 Distributia puterii disipate

Functionare magnetron dreapta

Detaliu

Fig.III.15 Puterea disipata in concentratorul de camp si in sarcina

Functionare magnetron sus si dreapta

CAPITOLUL IV Modelarea Fenomenelor de Camp Termic

Cuprins

IV.1 Introducere

IV.2 Determinarea profilelor de camp termic

IV.3 Programul informatic MacGFEM

IV.4 Concluzii

IV.1 Introducere

Modelarea fenomenelor de camp termic in cazul corpurilor procesate termic cu ajutorul microundelor, se face pornind de la cunoasterea precisa a campului electromagnetic in volumul acestuia.

Cunoscind acum valorile de camp , asa cum am vazut in CAP II, vom putea calcula puterea disipata in nodurile retelei-mash,

Desigur ca metodele variationale si metodele integrate permit rezolvarea problemelor de limite si discontinuitati, cu toate acestea metoda elementului finit este de multe ori preferata datorita rezultatelor bune obtinute in rezolvarea ecuatiilor de camp si a celor termice. Un alt avantaj al metodei elementului finit este ca poate fi folosita cu succes in cazul incintelor cu sarcini plasate aleator si forme geometrice complexe.

IV.2 Determinarea profilelor de camp termic

În ceea ce privește ecuația căldurii, care exprimă transferul de energie termică, prin conducție în interiorul unui volum încălzit cu microunde, aceasta se scrie:

Prin încălzirea cu microunde, chiar la echilibru termic, temperatura nu este constantă, datorită condițiilor la limită (sarcina nu este izolată termic) și datorită neomogenității puterii microundelor.

Pe suprafața corpului de probă, condițiile de limită sun aceleași cu cele din încălzirea clasică. În cazul în care corpul nu este izolat termic, există pe suprafață pierderi de energie prin conducție sau radiație în infraroșu. În cazul unei încălziri mixte, pe această suprafață va exista un câștig de căldură.

Lucrând în 2D, vom avea (x,y)+/t, luând intervale de timp mici, t=0,1ms, obținem pentru temperatură , o variație liniară.

Conform metodei elementului finit, volumul V, este segmentat într-o multitudine de domenii, definite de rețeaua de triunghiuri izoparametrice de ordinul I.

Condițiile la limită pe frontiera , sunt :

Expresiile lui cp, și k, se pun sub forma unor funcții, care să exprime valorile acestora în oricare regiune definită a structurii, astfel:

așa cum se poate vedea în programul utilizat .

Suntem în măsură să obținem acum, profilele termice în sarcină și în transformatorul de căldură.

Iată deci o prezentare comparativă a profilelor termice pentru cavitatea rezonantă în Fig.IV.1 până la Fig.IV.5, pentru următoarele situații :

l.Cavitate cu sarcină, fără transformator de căldură

2.Cavitate rezonantă, fără sarcină , cu transformator de căldură

3.Cavitate rezonantă, cu transformator de căldură și cu sarcină

– cazul cu 2 magnetroane

Fig IV.5 Profil termic- cavitate cu concentrator de camp si sarcina

Functionare magnetron sus si dreapta – t =33 min

O observație foarte importantă este de aceea că sarcina ceramică aflată în interiorul transformatorului de căldură, suferă o încălzire mixtă, datorită câmpului de microunde și prin radiație și convecție din partea transformatorului de căldură.

Este evident că datorită diferenței proprietăților absorbante față de microunde a ceramicii din transformatorul de căldură și a ceramicii industriale, ca și datorită ecranării față de câmp la care este supusă ceramica industrială, concentrarea câmpului electric este net superioară în transformatorul de căldură. Urmarea imediată va fi o creștere în temperatură diferențiată, accentuată și de proprietățile de conductivitate termică mai bune ale transformatorului de căldură, astfel încât acesta din urmă va urca în temperatură foarte rapid, apărând imediat un gradient de temperatură între cele două ceramici, ceea ce va duce la un transfer termic prin convecție și radiație între acestea.

Se poate aprecia că încălzirea ceramicii industriale se face preponderent prin convecție și radiație, față de încălzirea datorată microundelor.

Se poate observa că, cuplajul între temperatură și câmp, este chiar mai complex decât acela definit de sistemul de ecuații . Într-adevăr, evoluția câmpului și a temperaturii sunt legate de aplicator, iar simulările care se pot face deși nu dau rezolvarea exactă a problemei , ne oferă informații foarte interesante pentru înțelegerea încălzirii și pentru urmărirea întregului proces.

IV.3 Programul informatic MacGFEM

IV.4 Concluzii

Interpretind datele obținute se pot desprinde mai multe concluzii, astfel:

Plasarea sarcinii în cavitate, fără transformatorul de căldură nu dă rezultate satisfăcătoare din punct de vedere al intervalului de timp necesar pentru a ajunge la o temperatură dorită. Explicația ține de faptul că ceramicile au în general slabe proprietăți absorbante de microunde la temperatura camerei. De aceea, procesul de încălzire de pînă la l50C este lent.

Prin urmare soluția de obținere a unor temperaturi înalte, peste 11ooC, în timpi de ordinul minutelor, este aceea de plasare a sarcinii într-un creuzet concentrator de câmp avânnd proprietăți absorbante de microunde foarte puternice, chiar la temperatura normală. Acest obiect numit transformator de căldură este un material ceramic compozit, care are rolul de sursă de căldură, transmițând corpului de probă, energie termică prin radiații în infraroșu și prin convecție.[64]

Efectul este unul benefic pentru scopul urmărit, sarcina fiind supusă unei încălziri mixte de volum și de suprafață în același timp, cu rezultate bune privind calitatea produsului de sinterizare și ardere.

În ceea ce privește folosirea unuia sau a două generatoare de microunde, s-a putut constata că din punct de vedere al concentrării și repartizării câmpului, al puterii disipate și al temperaturilor este de preferat utilizarea simultană a două surse.

Deasemenea puterea nominală a surselor este optim aleasă, pentru valoarea de l kW, ținând cont de geometria structurii și de considerente economice.

Motivul alegerii metodei elementului finit în 2D și nu în 3D, își află explicația în dimensiunile mari ale cavității rezonante cu dificultățile majore ce apar privind generarea rețelei și a puterii de calcul necesare, de asemenea în problemele apărute din cauza folosirii transformatorului de căldură.

Cu ajutorul modelării în 2D am obținut rezultate pe care le putem aprecia bune, mai ales în ceea ce privește câmpul electromagnetic, rezultate care analizate comparativ au întărit ideea necesității utilizării acestuia și a folosirii simultane a două generatoare de microunde.

Modelarea poate oferi solutii rapide si putin costisitoare, atunci cind se incearca metode noi de lucru, sau cind se folosesc materiale cu un comportament necunoscut.

CAPITOLUL V Software utilizat în cercetarea din domeniul microundelor

Cuprins

V.1 Instrumente pentru modelarea și simularea câmpului magnetic, electric și termic

V.2 Analiza statistică a datelor experimentale

V.3 Software pentru studiul meta-analitic al literaturii de cercetare

V.4 Resurse disponibile în domeniu.

Softul utilizat în cercetarea din domeniul microundelor este dezvoltat pe doua direcții:

A. Instrumente pentru modelarea și simularea câmpului magnetic, electric și termic; aceste instrumente folosesc calculul numeric pentru rezolvarea ecuațiilor de câmp

B. Analiza statistică a datelor experimentale; ea permite elaborarea și testarea de ipoteze statistice.

În plus se mai utilizează și categoria:

C. Software pentru studiul meta-analitic al literaturii de cercetare, care sunt metode de căutare a publicațiilor relevante pentru o temă și analiza relevanței.

V.1 Instrumente pentru modelarea câmpului magnetic, electric și termic

Modelarea câmpului magnetic urmărește în principal 2 aspecte: analiza câmpului și simularea.

V.1.1 Clasificare

Softul de modelare a câmpurilor poate fi clasificat după mai multe criterii:

1. din punct de vedere matematic, după metoda de rezolvare a ecuațiilor cu derivate parțiale, în:

– FEM (Finite Element Method)/ FEA (Finite Element Analysis) utilizează Metoda elementului finit, care se bazează pe aproximări pe subdomenii, nefiind condiționată de forma domeniului (în spațiu)

– FDM (Finite Difference Method), Metoda diferențelor finite, care este mai lentă decât FEM [50, Zhao H. et al., 1998] și care se pretează cel mai bine pentru domenii rectangulare

-FDTD (Finite Difference Time Domain), Metoda diferențelor finite în domeniul de timp

-BEM (Boundary Element Method), Metoda elementelor de frontieră

-FVM (Finite Volume Method), Metoda volumelor finite

-FVTD (Finite Volume Time Domain); Metoda volumelor finite în domeniul de timp; o aplicație de uscare a lemnului cu microunde modelată pe FVTD poate fi consultată în [50, Zhao H. et al., 1998]

-VIM (Volume Integral Method), Metoda integralei de volum.

2. dupa accesibilitatea comercială: putem împărți produsele în soft gratuit și proprietar (proprietatea unei firme care vinde dreptul de utilzare);

3. dupa accesibilitatea în folosire: versiuni desktop, respectiv existența de module online și/sau cloud; versiunea desktop presupune instalarea programului în calculatorul utilizatorului, versiunea online presupupune că soft-ul se găsește pe site-ul proprietar, utilizatorul îl folosește de la distanță în mod individual, iar versiunea cloud presupune că mai mulți utilizatori pot colabora la simularea aceluiași model prin Internet.

În selecția pachetelor software din acest capitol am utilizat căutarea manuală combinată cu tehnica Ancestry Approach sau Descendancy Approach, care va fi prezentată în secțiunea 'C. Software pentru studiul meta-analitic al literaturii de cercetare', metodă utilizată și în [41, Popovici D., Popovici O., 2012];

Criteriile de evaluare utilizate pentru selecție au fost:

– tehnica numerică utilizată de soft

– disponibilitatea comercială (proprietar sau gratuit)

– disponibilitatea de folosire: dacă există module online.

V.1.2 FEM Comerciale/Proprietare

Cele mai cunoscute pachete software bazate pe Metoda elementului finit, care este și metoda cea mai frecvent utilizată la realizarea de software de modelare/ simulare, sunt:

MatLab

SimuLink

Flux 2D/3D

Comsol Multiphysics

AnSys Multiphysics

Concerto

QuickWave

Agilent EEsof EDA

Sonnet

MathCad

LabView

Ele sunt prezentate în continuare.

MATLAB® (MathWorks Inc., USA)

Denumirea produsului vine de la Matrix Laboratory. MatLab este un mediu de dezvoltare interactiv pentru calcule și grafică din domeniul tehnic. El folosește un limbaj de programare cu același nume și a fost creat de MathWorks Inc. [30, MathWorks, MatLab].

MatLab este un instrument de calcule numerice orientat spre matrici și vectori care integrează vizualizarea cu programarea. El are o interfață grafică intuitivă și poate interacționa cu alte aplicații prin import/ export de date (CAD-uri) sau prin interfețe (cu limbajele de programare C/C++, Fortran). Deasemenea utilizatorul poate să-și creeze propriile biblioteci de funcții, dezvoltând mediul.

MatLab este utilizat extensiv în industrie și universități și are variante pentru aproape toate platformele: Windows, GNU/Linux, UNIX și Mac OS.

Produsul are 5 componente de bază:

a) Mediul de lucru (GUI)

b) Biblioteca de funcții matematice

c) Handle Graphics, sistemul de grafică

d) Limbajul (de nivel înalt), dispune de instrucții procedurale și OOP

e) Interfața de aplicații program (API) permite scrierea de programe (C/C++, Fortran).

În afara modulului de bază, dispune de extensii orientate pe domenii specifice, identificabile sub denumirea de ToolBoxes, ca de exemplu:

– PDE (Partial Differential Equations) ToolBox, rezolvă ecuații diferențiale parțiale cu metoda elementului finit [28, MathWorks, 1984-1997];

– Signal Processing ToolBox, se folosește pentru analiză și reprezentare semnal;

– Statistic ToolBox conține: Distribuții de probabilități: Chi-square, Distribuții Poisson, Regresii, Analiza Multivariată, Analiza exploratorie a datelor, Testarea ipotezelor statistice [22, Ghinea M., Firețeanu V., 1998]

– Model Predictive Control ToolBox, rezolvă ecuații pentru controlul după un model dinamic explicit (o traiectorie), cum este cel al mișcării (robotică).

SimuLink® (MathWorks Inc., USA)

Un pachet aditional MatLabului și destinat simulării, este produs tot de firma MathWorks și distribuit sub denumirea de SimuLink. El oferă oportunitatea de a realiza modelarea și simularea sistemelor dinamice (mecanice, electrice) liniare și neliniare, continue și discrete, ca și a celor multi variabile, folosind modele matematice. Modelul matematic al unui sistem dinamic constă în reprezentarea dependențelor dintre componente. Dacă modelul corespunde unui proces fizic atunci este un model sistematic, pentru care există relații cauzale între cantități.

SimuLink este un simulator care are o interfață grafică în care modelele sunt prezentate ca blocuri și pot fi manipulate prin tehnica drag_and_drop, ceea ce simplifică mult programarea [23, Iliescu S., 2005]. SimuLink este convenabil mai ales pentru procesare semnalelor, ca și pentru proiectarea controllerelor ca PID (Proportional-Integral-Derivative) sau a controllerelor MPC (Model Predictive Controllers) utilizate în robotică.

Pe site-ul producătorului se găsește o bogată documentație care descrie toate modulele, atât ale MatLabului cât și ale SimuLinkului [29, MathWorks, Documentation Center].

Flux 2D/3D® (Cedrat, Franța)

Flux este o aplicație bazată pe Metoda elementului finit, având o interfață prietenoasă și oferind simulări electromagnetice și termice, ca și proiectare și optimizare. Este potrivit pentru analiză statică, armonică și a regimurilor tranzitorii ale câmpurilor magnetice, electrice și termice; analiză parametrizată, cuplarea câmpurilor magnetic/ electric/ termic, conectarea circuitelor externe, cuplare mecanică, în 2 și 3D. Dispune și de un modul Flux/ SimuLink [16, Cedrat, 1998].

Aplicații: mașini rotative, actuatoare lineare, dispozitive de încălzire prin inducție, transformatoare și inductanțe, senzori, dispozitive de înaltă tensiune, cabluri, compatibilitate electromagnetică, testare non distructivă.

Flux acceptă o descriere geometrică directă sau specificații CAD ca date de intrare, și are capabilități mesh și analiză, utilizând pentru discretizarea domeniului de calcul 2D o rețea triunghiulară cu 6 noduri (în mijloacele laturilor și în vârfuri), sau o rețea dreptunghiulară; pentru 3D poate combina elemente tetraedrale, piramidale și paralelipipedice; deține un solver rapid și robust, capabilități de post-procesare și export în CAD; rezultatele se obțin sub formă de valori locale și globale ale mărimilor, curbe de evoluție și hărți de culoare; un Trial este oferit de magsoft-flux.com [27, Magsoft Corporation].

Poate fi folosit în Modelarea numerică a câmpului electomagnetic în instalațiile de încălzire prin inducție sau în Procesarea materialelor dielectrice în câmp de radiofrecvență (modulul Dielectro-Thermal) pentru determinarea distribuției câmpului electromagnetic respectiv termic în diverse dispozitive, inclusiv dielectrici imperfecți care suportă tensiuni sinusoidale [49, Stoichescu, 2010];

Comsol Multiphysics® (Comsol, Suedia)

Comsol (care s-a numit inițial FEMLAB) este un software de analiză și simulare pentru procesele fizice care pot fi descrise prin sisteme de ecuații diferențiale, bazat pe metoda elementului finit. Poate fi utilizat pentru o multitudine de fenomene și aplicații inginerești, mai ales pentru frecvențe mari și la fenomene de cuplare de câmp (mai puțin câmpul de masă).

Are un mare număr de module pentru programare, preprocesare și postprocesare, precum și o interfață pentru MatLab (Live Link). Include modul pentru ecuații cu diferențe parțiale ce pot fi introduse direct sau în formă slabă, pentru modelarea performațelor condensatoarelor, pentru transferul de căldură (conducție, convecție, radiație), modul radiofrecvență, import geometrii CAD, modul scripting (pentru limbajul propriu care permite extinderea analizei, a reprezentărilor grafice și realizarea de interfețe personalizate) șa. Este disponibil pe toate platformele uzuale (Windows, Mac, Linux) [9, Comsol, 2008].

Poate fi folosit în Procesarea materialelor dielectrice în câmp de microunde, în Industria lemnului sau în prelucrarea alimentelor, alte procese industriale de uscare (ex. circuite integrate), mecanica fluidelor, probleme cuplate.

Firma oferă Work-shop-uri, traininguri, precum și Webinare gratuite, care sunt sesiuni online de prezentare și de training pentru produs, pe diferite teme, ca de exemplu: Simularea transferului de căldură, Modelarea transformatoarelor și inductoarelor, Simularea în Biotehnologie, Utilizarea MatLab cu Comsol, Inginerie Biomedicală, Încălzire cu Radiofrecvențe și cu Microunde șa. Webinarele sunt înregistrate și pot fi audiate și ulterior difuzării lor, la cerere [10, Comsol, 2013].

ANSYS Multiphysics (AnSys, SUA)

Este un produs puternic, având următoarele componente [24, INAS]:

– Structural Mechanics: simulare în domeniul structural/ termic, simulează comportarea diferitelor materiale;

– Fluid Dynamics: simularea curgerii fluidelor;

– Explicit Dynamics: studiul comportamentului dinamic neliniar al materialelor;

– Electromagnetics LF: studiul câmpului electromagnetic de joasă frecvență, de la comunicații mobile la până la sisteme electromecanice;

– HFSS High Frequency Structure Simulator (AnSys, SUA);

Componenta se adresează în special proiectării de antene și proiectării de circuite electronice complexe în radiofrecvență, ca și pentru echipamente biomedicale [2, AnSys HFSS]; are două moduri de operare: modul 3-D utilizat pentru modelare sau import CAD a geometriilor complexe și simulare de înaltă frecvență; cel de-al doilea mod, HFSS pentru ECAD, este pentru proiectanții care lucreză în geometrii pe straturi; HFSS are un solver bazat pe metoda elementului finit care se poate folosi în domeniul de înaltă frecvență, regim tranzitoriu, ecuații integrale și optică.

HFSS a fost dezvoltat de compania AnSoft, care a fost ulterior cumpărată de AnSys [1, AnSys HFSS].

O comparație între Comsol și AnSys a publicat [47, Salvi D. et all, 2010], pentru simularea încălzirii cu microunde a unui lichid în curgere continuă, concluzionând că rezultatele obținute cu cele două pachete sunt comparabile; ca diferențe ei au remarcat că stabilirea modelului pare mai flexibilă cu Comsol, dar în același timp acest produs consumă mai multă memorie.

Concerto (AVL, Austria)

Producătorul soft-ului, compania AVL, s-a ocupat cu precădere de motoarele diesel, aducând numeroase contribuții la dezvoltarea acestora. Modelatorul Concerto este disponibil în mai multe pachete predefinite în funcție de necesitățile utilizatorilor [11, Concerto AVL]:

– Cocerto- Key este pachetul de bază și conține post-procesare generică (de ex. diagrame 2D, obiecte tabelă, obiecte text, funcții formule, acces la formate de fișiere de bază);

– Concerto-Advanced for Indicating, include funcționalitatea completă a Concerto-Key precum și funcții suplimentare relevante pentru sarcini de analiza combustiei (ex. editorul de formule grafice CalcGraf cu biblioteca de formule pentru analiza combustiei);

– Concerto-Advanced for Test Bed, include funcționalitatea completă a Concerto-Key precum și funcții suplimentare relevante pentru evaluarea test bed (ex. acces la sistemele de baze de date ASAM-ODS, suport pentru engine de hărți 3D);

– Concerto-Premium, include funcționalitatea completă a Concerto-Key, Concerto-Advanced for Indicating, Concerto-Advanced for Test Bed și adițional un set de funcții expert (ex. access extins la date, interfețe COM- și Active-X pentru control avansat de la distanță, funcții avansate de căutare a datelor, adaptabilitatea interfeței cu utilizatorul) [12, Concerto, 2010].

QuickWave-3D (QWED, Polonia)

Este un solver/ simulator FDTD (Finite Difference Time Domain), destinat proiectării electromagnetice și măsurării proprietăților electromagnetice ale materialelor, având și o selecție de modele unice pentru frontiere curbe, interfețe de mediu, excitație modală și extragerea parametrilor [43, QuickWave].

QW conține două programe principale:

– QW-Editor, pentru definirea geometriei, a mediului, parametrii de Intrare/iEșire și postprocesare; mai conține o bibliotecă cu obiecte parametrizate și capacitatea de a genera alte obiecte și biblioteci; conversia la/ de la formatul CAD este facilitat;

-QW-Simulator, este un solver unic, FDTD; output-urile lui includ matrici S multi-modale și multi-port; norul de puncte al radiației, aspectul câmpului, puterea dispată, reflectometria pe domeniul timp șa.

Printre modulele opționale se numără:

– QW-BHM (QuickWave- Basic Heating Module) specializat pe încălzirea cu microunde, ce include analiza încălzirii corpurilor în mișcare: QW-HFM (QuickWave Heat Flow Module)

– QW-AddIn for Autodesk Inventor Software, plug-în-ul pentru Inventor, ce permite definirea unei simulări electromagnetice complete direct în Autodesck Inventor utilizând convenabila interfață QW

– QW-GPUSim, un limbaj OpenCL folosit pentru implementarea Simulatorului QW, care rulează pe plăcile grafice moderne și care permite accelerarea calculelor aproximativ de 10 ori la modelele 3D și de 20 de ori la modelel V2D, păstrând multithreading-ul

– S-Converter, modul specializat pe calcularea parametrilor întregii structuri pornind de la parametrii unei părți, măsurați sau calculați la un port extern al structurii.

ADS Agilent EEsof EDA (Agilent Technologies, USA)

Sub denumirea de AppCAD se cunosc 2 pachete software:

AppCAD de la Avago Technologies, care este un produs Free și este prezentat la secțiunea respectivă, în timp ce grupul de aplicații 'AppCAD' de la Agilent Technologies este un grup de produse soft orientate mai ales pe analiza de semnal, dezvoltat de o companie care este și producătoare de echipamente: osciloscoape, surse șa. [3, Agilent]

Agilent Technologies este fosta divizie de teste, măsurători, semiconductoare și analiză chimică de la Hewlett-Packard și a preluat AppCAD cu ocazia divizării HP, ulterior secțiunea de semiconductoare s-a separat încă odată în Avago Technologies, preluând AppCAD [4, AppCAD].

Pachetul de aplicații proprietare dezvoltat ulterior de Agilent conține 17 aplicații cu MatLa ca bonus, iar ADS (Advanced Design System) este aplicația orientată pe proiectarea circuitelor în RF. Agilent EEsof EDA software permite proiectarea, modelarea și simularea echipamentelor de înaltă frecvență.

SONNET (Sonnet Software Inc., USA)

Este o suită comercială de analiza circuitelor uni și multistrat planare, și a antenelor. Folosește Metoda Momentelor aplicată direct Ecuațiilor lui Maxwell pentru rezolvarea problemelor planare [48 SONNET]. Pe site-ul producătorului se găsește o descriere amănunțită a modului în care softul lucrează.

O versiune limitată, Sonnet Lite, este oferită gratuit.

PCT MathCad (PCT, USA)

Se adresează calculelor inginerești din domeniul Construcții civile, Inginerie mecanică și Inginerie electrică, incluzând capabilitatea de transformare a unităților de măsură [39, PCT MathCad]; este conceput sub formă de 3 biblioteci interactive pe domeniile de mai sus și 4 biblioteci extensii:

– MathCad Data Analysis Extension Pack

– MatcCad Signal Processing Extension Pack

– MathCad Image Processing Extension Pack

– MathCad Wavelets Extension Pack;

Pe MathCad Academic Channel (YouTube) se găsește o bogată colecție de tutoriale video destinată învățării soft-ului.

LabView (Laboratory Virtual Instrumentation Engineering Workbench National Instruments, USA)

O aplicație pentru instrumentație virtuala este LabView. Se folosește pentru creerea de echipamente de măsură și control simulate, inclusiv interfețe, realizându-se astfel economii considerabile;

LabView în România: http://romania.ni.com/labview

V.1.3 FEM Free

Cele mai cunoscute produse Free, și care vor fi prezentate mai jos sunt:

FreeFEM

Mefisto 2D/3D

AppCAD

Aurora Z88

OCCT& Salome

FreeFEM (licență GNU/GPL, Franța)

FreeFEM este un solver de ecuații cu diferențe parțiale (EDF) prin metoda elementului finit, care include un limbaj propriu de programare de nivel înalt. Rezolvă EDF în spațiul 2/3D, folosind o mare varietate de elemente finite, inclusiv elemente discontinue. Este cel mai puternic solver FEM distribuit sub licență GPL, ceea ce înseamnă că utilizatorul poate avea acces și la sursa programului, pe care poate s-o modifice, cu condiția să o distribuie mai departe [18, FreeFEM].

Istoricul produsului poate fi împărțită în 3 etape: Inițial a existat o varintă pentru MacIntosh sub denumirea de MacFem, utilizată și de Popovici O. [42, Popovici O., Popovici D., 2008] la modelarea câmpului de microunde; Între 1992-2012 au fost puse la dispoziția utilizatorilor 4 variante: FreeFEM (1992), FreeFEM+ (1996), FreeFEM3D și FreeFEM++ v3.2 [19, FreeFEM]; din septembrie 2013 este disponibilă versiunea FreeFEM++ v 3.25 pentru Mac OS X și v 3.23 pentru Windows.

Există un proiect pentru utilizarea online a FreeFEM++ și a FreeFEM3D, prin intermediul interfeței FFW (FreeFEM on the Web), conform [20, FreeFEM], dar aceasta încă nu funcționează.

Grupul dezvoltatorilor are un nucleu la Universitatea Pierre et Marie Curie, Paris, Franța, în jurul profesorului Frédéric Hecht. La dezvoltarea acestui produs a contribuit și românul Ionuț Dănăilă, doctor în mecanica fluidelor, actualmente profesor la Universitatea din Rouen, Franța [15, Dănăilă I.].

Ionuț Dănăilă este autor, alături de Frédéric Hecht și Olivier Pironneau (UPMC) a cărții Simularea numerică în C++, Ed. Dunod, 2003, care conține exemple de utilizare a FreeFEM. În 2005 Ionuț Dănăilă și colaboratorii, au scos încă o carte de calcul științific în MatLab (tot la editura Dunod), tradusă în limba engleză și publicată de editura Springer în 2007; cărțile sunt prezentate în secțiunea D. Resurse.

El a luat parte la Work-shop-urile anuale 'FreeFem workshop on Generic Solver for PDEs: FreeFem++ and Applications', din perioada 2009- 2012, Paris, Jussieu [21, FreeFEM], care au avut obiective ca: Condensatele Bose-Einstein sau Gradienții Sobolev. În prezent lucrează la proiectul Becasim (Bose-Einstein Condensates: Advanced SIMulation deterministic and stochastic computational models) [8, Becasim], care are ca obiectiv dezvoltarea unui solver bazat pe calculul paralel.

Mefisto 2D/3D (licență GNU/GPL , Franța)

Mefisto 2D/3D este un produs OpenSource realizat tot în grupul 'Jussieu', de la Universitatea Pierre et Marie Curie, Paris, de Alain Perronnet împreună cu Pascal Joly și Christophe Doursat. El poate rezolva Ecuații Diferențiale Parțiale prin Metoda Elementului Finit, și este scris în Fortran 77 și C. Se poate folosi pentru simularea propagării microundelor ca și pentru proiectarea, analiza și studiul circuitelor cu microunde.

Are următoarele componente [31, Mefisto 2D/3D, Content]:

– Mefisto-Initier și Mefisto-Mesher permit creerea de mesh 2D/3D pornind de la puncte, linii, suprafețe, volume; construcția unui obiect cu interpolarea Lagrange isoparametrică a unor polinoame de gradul 1 și 2, respectiv cu metodologia de construcție Mefisto;

– Mefisto-Heater permite calculul ecuațiilor de transfer al căldurii și ecuațiile de undă de gradul 2; sunt disponibile încălzirea axi-simetrică și solver de elasticitate

– Mefisto-Elasticier calculeaza ecuațiile termo-elastice

– Mefisto-Fluider calculează ecuațiile Navier-Stokes constante și variabile.

Mefisto poate fi descărcat de la [32, Mefisto 2D/3D, Download].

AppCAD Free Design Assistant Tool (Avago Technologies, USA)

Avago este o companie americană specializată pe producție, în special pe mouse optic, dar și în domeniul de radio frecvență și microunde [5, AppCAD, Avago]; AppCAD free se adresează mai ales proiectanților de dispozitive cu circuite integrate și are următoarele module:

– Analiza și plotarea S-parametrului

– Proiectarea Circuitelor Active

– Analiza zgomotului în cascadă și a IP3

– Analiza liniilor de transmisie

– Semnale și Sisteme

– Modul de calcule inginerești complexe.

O video prezentare, filmată la standul Avago Technologies de la International Microwave Symposium IMS 2012, Canada, poate fi văzută la [6, AppCA, Avago].

Aurora Z88 (licență GNU/GPL, Germania)

Este un bun FEM solver promovat de Universitatea din Bayreuth, Germania. Este destinat calculelor statice lineare în mecanică, pe baza Metodei elementului finit. Poate fi utilizat și analiza de termie statică și în analiza frecvențelor naturale [7, Aurora Z88].

Open CASCADE Technology OCCT& Salome (OpenCasCade SAS, Franța)

Este o platformă OpenSource de dezvoltare în C++ și Pyton, pentru modelare suprafețe și solide, vizualizare și simulare [34, OCCT], dezvoltată de un grup francez din Guyancourt în colaborare cu cercetători ruși din Nizhny Novgorod; la origine a aparținut grupului Matra Datavision, preluat de IBM; în prezent OpenCasCade SAS este o subsidiară a Areva Group.

OCCT servește în principal la modelare iar Salome [46, Salome] se adresează utilizatorilor de pre și post-procesare pentru softul de analiză cu elemente finite.

O listă cu alte solvere PDE Free poate fi găsită la [36, OpenNovation], respectiv o listă cu alte soft-uri FEM, free, opensource sau proprietare, poate fi găsită la [52, wikipedia].

V.2 Analiza statistica a datelor experimentale

V.2.1 Utilizare software de statistică

Software-ul de statistică se potrivește fenomenelor probabilistice, cum ar fi de exemplu estimarea randamentului unei instalații cu microunde. Industria software produce un mare număr de pachete, mult mai mare ca în domeniul modelării și simulării de câmp, astfel încât înainte de a alege un produs avem nevoie de o documentare temeinică.

A apărut astfel necesitatea comparațiilor între performanțele acestor pachete. În alegerea pachetelor discutate ne-am orientat și după alte sudii publicate [40, Popovici D., 2013] și [51, Zhu X., 2005].

V.2.2 Pachete comerciale/proprietare

Origin

Origin (OriginLab, SUA)

Cel mai folosit pachet statistic în SUA, are 3 componente:

– analiza de date (regresii liniare, polinomiale, multiple și alte curbe, procesare de semnal și imagine, calcul științific, statistică descriptivă, teste non și parametrice,

– reprezentare grafică diferite forme

– programare cu LabTalk, limbaj de scripting nativ, cu Origin C limbaj compatibil ANSI C, X-Functions permite creerea de rutine proprii, Automation Server permite aplicații client LabView, Excel, MatLab, sau rutine scrise în Visual Basic, Visual C++, .NET șa [37, Origin Lab].

Pachetul se compune din Origin, OriginPro care conține extensii profesionale în statistică, procesare de semnal și imagine, și varianta cu preț redus (începe de la 70 USD) OriginPro [anonimizat] 1-3 ani care se distribuie studenților etc.. Tutoriale se pot găsi și pe YouTube, ca și [38, Origin].

V.2.3 Pachete Free

R-Project

Instat

Un studiu pentru pachete free mai ușoare, care se pot utiliza în scop educațional a publicat [51, Zhu X.& Kuljaca O., 2005]; ei au comparat următoarele suite: Visual Statistics ViSta (prof. Forrest Young, University of North Carolina, USA), AM (American Institute for Research, USA), WinIDAMS (UNESCO), IRRISTAT (International Rice Research Institute, Filipine), OpenStat3 (former prof. William Miller, Iowa State University, USA), PAST (Natural History Museum, Norvegia) și Instat+ (University of Reading, Marea Britanie); Criteriile de selecție folosite de ei, au fost:

– pachetele să fie gratuite

– să fie ușor de învățat și ușor de utilizat, astfel încât studenții să se poată concentra pe procesarea datelor și nu pe învățarea instrumentului

-să aibă tutoriale bune și componente pentru învățare (teaching components)

– să conțină testele comune cum ar fi: Chi pătrat, tabele de frecvență, analiză univariată și multivariată, ANOVA- analiză de varianță, box-plots, diagrama nor de puncte și regresie liniară

– să fie compatibile Windows.

R-Project (licență GNU/GPL, Auckland University, Noua Zeelandă)

R-Project este un proiect GNU. GNU a fost inițial o încercare de a produce o clonă Unix gratuită, din care a rezultat sistemul de operare Linux; denumirea 'GNU' este un acronim recursiv pentru 'GNU's Not Unix', ceea ce vrea să sublinieze că noul nucleu este diferit de nucleul Unix-ului; promotorul proiectului, Free Software Foundation, a elaborat ulterior specificațiile General Public License (GNU/GPL) prin care a stabilit principiile/ conceptul de software liber de a fi utilizat.

Astfel, utilizatorii care primesc software GNU au datoria să-l distribuie mai departe în aceleași condiții în care l-au primit, ceea ce implică în esență interzicerea de a comercializa un software care a fost primit gratuit, chiar și modificat.

R-Project sprijină o serie de Conferințe organizate anual de membrii comunității useR! care au ca obiect atât probleme specifice R cât și probleme de statistică în general, cu sesiuni care se desfășoară în format liber și nu sunt urmate de publicarea Proceendings-urilor.

De asemenea R-Project sprijină și o platformă pentru dezvoltatori, DSC Directions in Statistical Computing, care a organizat și ea o serie de Conferințe bi-anuale, între 1999- 2009, al căror rezultat se publică în revista Computational Statistics (Springer Link ). Ambele conferințe sunt coordonate de R Foundation Conference Committee, RFCC.

O altă revistă, Journal of Statistical Software, redactată de un grup de voluntari de la UCLA Statistics, publică o serie de volume legate de R, ca de exemplu: Econometrie în R, Ecologie și Modelare ecologică în R, Imagerie cu Rezonanță magnetică în R sau Spectroscopie și Chemometrie în R.

R este un pachet 'greu', necesită un efort serios de învățare, în schimb ca și calitate, a fost de exemplu comparat sub aspectul exactității, al încrederii în rezultate, de [26, Keeling K., Pavur R., 2007], pentru versiunea 1.9.1, și de [35, Odeh O. et all, 2010], pentru versiunea 2.20, obținând rezultate foarte bune. În aprilie 2013 a apărut versiunea 3.0.0 a produsului [44, R-Project].

În ceea ce privește popularitatea acestui produs, Muenchen Robert de la Statistical Consulting Center, University of Tennessee, USA, a publicat pe site-ul personal [33, Muenchen R., 2013], date potrivit cărora:

– într-un sondaj al Rexer Analytics din SUA, cu întrebarea 'Ce instrument de Data Mining/ Analiză ați utilizat în 2009?', R conducea urmat de SAS și IBM SPSS

– R se găsea în 2011 pe primul loc între pachetele utilizate în competiții de Analiză a datelor, urmat de MatLab și SAS

– la numărul utilizărilor în publicațiile academice, măsurat prin Google Scholar, conducea SPSS urmat de SAS.

La: [45, R-Project, Tutoriale] pot fi găsite cîteva tutoriale împreună cu ebook-ul ' R Tutorial with Bayesian Statistics Using OpenBUGS', de Chi Yau. De asemenea un site plin cu informații utile începătorilor și utilizatorilor a fost creat de Rob Kabacoff, autorul cărții 'R in action: Dana analysis and Graphics with R', publicată în 2013 la a doua ediție.

Instat (University of Reading, Marea Britanie)

Instat dispune de interfață grafică GUI, și deși are unele capacități de plotare discutabile, este foarte bun pentru predare/ învățare, aceasta fiind de altfel una dintre intențiile autorilor. Toate funcțiile uzuale sunt prezente, incluzand analiza multivariată, dar fără regresia nonlineară. De altfel, majoritatea pachetelor gratuite nu suportă tehnicile nonlineare.

Instat este gratuit pentru aplicații noncomerciale monoutilizator, dar necesită costul unei licențe pentru alte tipuri de utilizatori [25, Instat]. Se bazează foarte mult pe asemănarea cu foile de calcul tabelar, dar nu permite programarea (cu excepția macrourilor). Instat+ este versiunea pentru Windows.

[51, Zhu X.& Kuljaca O., 2005] afirmă că este depășit ca funcționalitate de ViSta și WinIdams.

Alte pachete sunt descrise în [40, Popovici D., 2013]: SPSS (IBM, USA), Statistica (StatSoft Company, USA), Minitab (Minitab Inc., USA), GraphPad (GraphPad Software Inc., USA), Gpower (University of Dusseldorf, Germania) și Epi Info (Center for Disease Control, Georgia, SUA).

V.3 Software pentru studiul meta-analitic al literaturii de cercetare

Volumul literaturii de specialitate din toate domeniile a crescut în ultimul timp foarte mult, de aceea sunt necesare tehnici de căutare a publicațiilor relevante în Bazele de date. Tehnica Ancestry Approach sau Descendancy Approach folosește referințele bibliografice are articolelor/ cărților identificate, mergând din aproape în aproape [13, Cooper H., 1982].

O atenție sporită merită acordată tezelor de doctorat, ca fiind purtătoare de elemente de noutate/ actualitate.

V.4 Resurse disponibile în domeniu

International Microwave Power Institute IMPI, USA (http://www.impi.org/Home_Page.html) publica o prestigioasă revistă trimestrială: The Journal of Microwave Power and Electromagnetic Energy (JMPEE) (http://www.jmpee.org/), indexată de Thomson Reuters în Science Citation Index Expanded (SciSearch), Journal Citation Reports/Science Edition și Current Contents/Engineering Computing and Technology. Revista poate fi consultată și pe portalul Research Gate: https://www.researchgate.net/journal/0832-7823_The_Journal_of_microwave_power_and_electromagnetic_energy_a_publication_of_the_International_Microwave_Power_Institute

International Microwave Symposium, este organizat anual, pe continentul Nord American, de IEEE Microwave Theory and Techniques Society's, în cadrul Săpămânii Microundelor http://ims2012.mtt.org/

Association for Microwave Power în Europe for Research and Education (AMPERE, http://www.ampereeurope.org/), organizație non-profit cu sediul în Franța, organizează o Conferință anuală itinerantă, în fiecare an la altă universitate europeană. În 2007 Conferința, a 11-a ediție Ampere, a avut loc la Universitatea din Oradea, România.

Robert F. Schiffmann, State of the art of microwave applications în the food industry în the USA, Advances în Microwave and Radio Frequency Processing ed., M. Willert Porada (Ed.), Ed. Berlin, Germany: Springer-Verlag, 2006; Robert Schiffmann este membru în Comitetul executiv al IMPI

Dănăilă I., Joly P., Kaber S.M., Postel M., 2007, An Introduction to Scientific Computing. Twelve computational projects solved with Matlab, Springer

Dănăilă I., Joly P., Kaber S.M., Postel M., 2005, Introduction au calcul scientifique par la pratique : 12 projets résolus avec Matlab, Dunod (lb. franceză).

Bibliografie Capitolul V

Bibliografie

Anexe

Anexa A1 Moduri proprii unei structuri rectangulare

Rezolvarea ecuațiilor lui Maxwell, ținând cont de condițiile la limită impuse de pereții metalici, arată că într-un ghid de undă se pot propaga un număr infinit de diferite tipuri de unde electromagnetice, numite moduri. Aceste moduri formează o bază de descompunere a câmpului, într-un ghid de undă omogen, existând două tipuri de moduri:

1). Modurile TEmn

Modurile transvers electric sunt caracterizate prin prezența câmpului electric doar în secțiunea transversală pe direcția de propagare, câmpul magnetic având și o componentă longitudinală:

(A1-1)

(A1-2)

unde și , sunt constante de normare care se pot determina cu ajutorul relațiilor (A1-8) și (A1-9).

2). Modurile TMmn

Modurile transvers magnetic, sunt caracterizate prin prezența câmpului magnetic doar în secțiunea transversală pe direcția de propagare, câmpul electric având și o componentă longitudinală:

(A1-3)

(A1-4)

unde și , sunt constante de normare care se pot determina cu ajutorul relațiilor (A1-8) și (A1-9).

Ghidurile multimod, permit propagarea mai multor moduri de ordin superior, astfel, pentru ca un mod (m, n) să se propage, acesta trebuie să verifice condiția:

(A1-5)

Expresiile câmpurilor rămân aceleași ca pentru ghidul monomod.

Modurile transvers electric TE, sau transvers magnetic TM, constituie o bază completă de funcții.

Fie Γq, constanta de propagare a modului q:

(A1-6)

unde q = (m, n)

Notând cu și respectiv componentele transverse ale modurilor lui și a lui acestea verifică relația:

unde -sensul de propagare

-impedanța modului q

Pentru un mod TE:

Pentru un mod TM:

Condiția de propagare fiind:

Pentru a determina constantele de normare, vom aplica următoarea relație de ortogonalitate:

(A1-7)

Relația (A1-7), se poate scrie sub forma a două relații:

(A1-8)

(A1-9)

unde

Anexa A2 Fotografii echipamente de incalzire cu microunde

Cuptor cu microunde, cu sarcina ceramica, la 1050 grd C

Echipamente de automatizare-proces de incalzire

Aparate de masura-proces de incalzire

Procesarea informatica a datelor

Anexa A3 Simboluri

E – vectorul câmp electric

H – vectorul câmp magnetic

D – vectorul inducție electrică

B – vectorul inducție magnetică

J – densitate de curent

– conductivitate electrică

0 – permitivitate vidului

r – permitivitatea relativă a mediului

– densitatea de sarcină

’ – partea reală a permitivității complexe

” – partea imaginară a permitivității complexe

tg – tangenta unghiului de pierderi

– lungimea de undă

0 – lungimea de undă în spațiul liber

0 – constanta de atenuare

0 – constanta de fază

– permeabilitatea magnetică

k – număr de undă

– pulsația

f – frecvența

f0 – frecvența de rezonanță

P – puterea

Pi – puterea incidentă

Pr – puterea reflectată

Pabs – puterea absorbită

R – coeficient de reflexie

– impedanță de undă

Q – cantitate de căldură

m – masa

cp – căldura specifică

L – căldura latentă

h – constanta lui Plank

H – coeficient global de transfer termic

h – vector unitar normal

TE – transvers electric

TM – transvers magnetic

TEM – transvers electromagnetic

G – funcția lui Green diadică

im – constanta de propagare a modului “m”, în ghidul “i”

A – matricea A

– matrice de difracție

(*) – desemnează complexul conjugat

– randament

t – temperatura

c – viteza luminii

Du – densitatea de umplere

Nc – număr de coordonare

Z – parametru de măsurare a sfericității

r – tensiunea de suprafață a lichidului

– unghiul de contact

a – raza porului

dx – derivata in functie de x

W – umiditate

Rg – constanta gazului

Vm – volumul molar

pa – presiunea aplicată

plot – afiseaza

solve – rezolva

D – difuzivitate

C – concentrație

[D] – concentrație de defecte

Qo – coeficient de supratensiune, cavitate goală

Qc – coeficient de supratensiune, cavitate cu eșantion

c – impedanța de undă în spațiul liber

Rm – rezistența superficială a metalului

Dm – adâncimea de pătrundere a undei în metal

m,n,p – indici de mod

Bibliografie Capitolul V

Bibliografie

Anexe

Anexa A1 Moduri proprii unei structuri rectangulare

Rezolvarea ecuațiilor lui Maxwell, ținând cont de condițiile la limită impuse de pereții metalici, arată că într-un ghid de undă se pot propaga un număr infinit de diferite tipuri de unde electromagnetice, numite moduri. Aceste moduri formează o bază de descompunere a câmpului, într-un ghid de undă omogen, existând două tipuri de moduri:

1). Modurile TEmn

Modurile transvers electric sunt caracterizate prin prezența câmpului electric doar în secțiunea transversală pe direcția de propagare, câmpul magnetic având și o componentă longitudinală:

(A1-1)

(A1-2)

unde și , sunt constante de normare care se pot determina cu ajutorul relațiilor (A1-8) și (A1-9).

2). Modurile TMmn

Modurile transvers magnetic, sunt caracterizate prin prezența câmpului magnetic doar în secțiunea transversală pe direcția de propagare, câmpul electric având și o componentă longitudinală:

(A1-3)

(A1-4)

unde și , sunt constante de normare care se pot determina cu ajutorul relațiilor (A1-8) și (A1-9).

Ghidurile multimod, permit propagarea mai multor moduri de ordin superior, astfel, pentru ca un mod (m, n) să se propage, acesta trebuie să verifice condiția:

(A1-5)

Expresiile câmpurilor rămân aceleași ca pentru ghidul monomod.

Modurile transvers electric TE, sau transvers magnetic TM, constituie o bază completă de funcții.

Fie Γq, constanta de propagare a modului q:

(A1-6)

unde q = (m, n)

Notând cu și respectiv componentele transverse ale modurilor lui și a lui acestea verifică relația:

unde -sensul de propagare

-impedanța modului q

Pentru un mod TE:

Pentru un mod TM:

Condiția de propagare fiind:

Pentru a determina constantele de normare, vom aplica următoarea relație de ortogonalitate:

(A1-7)

Relația (A1-7), se poate scrie sub forma a două relații:

(A1-8)

(A1-9)

unde

Anexa A2 Fotografii echipamente de incalzire cu microunde

Cuptor cu microunde, cu sarcina ceramica, la 1050 grd C

Echipamente de automatizare-proces de incalzire

Aparate de masura-proces de incalzire

Procesarea informatica a datelor

Anexa A3 Simboluri

E – vectorul câmp electric

H – vectorul câmp magnetic

D – vectorul inducție electrică

B – vectorul inducție magnetică

J – densitate de curent

– conductivitate electrică

0 – permitivitate vidului

r – permitivitatea relativă a mediului

– densitatea de sarcină

’ – partea reală a permitivității complexe

” – partea imaginară a permitivității complexe

tg – tangenta unghiului de pierderi

– lungimea de undă

0 – lungimea de undă în spațiul liber

0 – constanta de atenuare

0 – constanta de fază

– permeabilitatea magnetică

k – număr de undă

– pulsația

f – frecvența

f0 – frecvența de rezonanță

P – puterea

Pi – puterea incidentă

Pr – puterea reflectată

Pabs – puterea absorbită

R – coeficient de reflexie

– impedanță de undă

Q – cantitate de căldură

m – masa

cp – căldura specifică

L – căldura latentă

h – constanta lui Plank

H – coeficient global de transfer termic

h – vector unitar normal

TE – transvers electric

TM – transvers magnetic

TEM – transvers electromagnetic

G – funcția lui Green diadică

im – constanta de propagare a modului “m”, în ghidul “i”

A – matricea A

– matrice de difracție

(*) – desemnează complexul conjugat

– randament

t – temperatura

c – viteza luminii

Du – densitatea de umplere

Nc – număr de coordonare

Z – parametru de măsurare a sfericității

r – tensiunea de suprafață a lichidului

– unghiul de contact

a – raza porului

dx – derivata in functie de x

W – umiditate

Rg – constanta gazului

Vm – volumul molar

pa – presiunea aplicată

plot – afiseaza

solve – rezolva

D – difuzivitate

C – concentrație

[D] – concentrație de defecte

Qo – coeficient de supratensiune, cavitate goală

Qc – coeficient de supratensiune, cavitate cu eșantion

c – impedanța de undă în spațiul liber

Rm – rezistența superficială a metalului

Dm – adâncimea de pătrundere a undei în metal

m,n,p – indici de mod

Similar Posts