SIMULINK

GEAVANCEERD SIMULEREN MET MATLAB / FEMLAB / SIMULINK ir. A.W.M. (Jos) van Schijndel, Technische Universiteit Eindhoven Faculteit Bouwkunde, Vakgroep FAGO Email: [email protected]

SAMENVATTING FemLab wordt geëvalueerd als gereedschap voor het oplossen van bouwfysische en installatietechnische problemen die gebaseerd zijn op partiele differentiaal vergelijkingen (PDEs). De software is ontworpen om systemen van gekoppelde PDEs, 1D, 2D of 3D, al dan niet lineair en tijdsafhankelijk te simuleren. Om te laten zien hoe het programma werkt, is een complete code alsmede de resultaten gegeven, voor het simuleren van een 2D stationair koudebrug probleem. Voor validatie doeleinden zijn simulaties van een 1D dynamisch niet lineair vochttransport probleem in een poreus materiaal en een 2D dynamisch luchtstromingsprobleem gebruikmakend van de Navier Stokes vergelijkingen en convectie, vergeleken met metingen. Alle simulatie resultaten komen goed overeen met de metingen. Om te laten zien hoe krachtig de FemLab software is, ook in combinatie met Simulink, worden de volgende applicaties behandeld: 1) het modeleren van 3D gecombineerd warmte- en vochttransport in een materiaal, 2) het modeleren van koudeval compensatie door een (put)convector, 3) het dynamisch simuleren van een 2D koudebrug in combinatie met een gebouwmodel door de 'standaard' export mogelijkheid van een FemLab model naar Simulink, 4) het regelen van een luchtstromingsmodel in FemLab met een aan/uit regelaar van Simulink. De FemLab software is geschreven in de Matlab simulatie omgeving en daardoor is het mogelijk gebruik te maken van de reeds aanwezige visualisatie tools, toolboxen (waaronder Simulink) en alle andere programma's geschreven in MatLab. Hierdoor is de combinatie MatLab/FemLab/Simulink een zeer krachtige en flexibele simulatie omgeving voor het oplossen van wetenschappelijke en technische problemen.

1. INTRODUCTIE Veel wetenschappelijke en technische problemen in de bouwfysica en installatietechniek kunnen beschreven worden met PDEs. Er zijn op dit moment veel software programma's beschikbaar waarbij een specifiek PDE probleem wordt opgelost. Deze software programma's zijn ontworpen om in korte tijd simulatie resultaten te verkrijgen en meestal ligt hierbij de nadruk op het eenvoudig kunnen invoeren van de data zoals bijvoorbeeld een geometrie en materiaaleigenschappen. Een nadeel is dat de software meestal niet flexibel is als de gebruiker een model wil veranderen of combineren. Hoe de gehanteerde modellen gesimuleerd worden en wat de kwaliteit van de numerieke oplossing is, is vaak niet bekend, daardoor is het rekengedeelte ook vaak een 'black-box'. Een andere categorie van software, zoals FemLab [1], FlexPDE [2] en PDEase2 [3], is speciaal ontwikkeld voor het simuleren van PDEs waarbij de gebruiker in principe ieder stelsel van gekoppelde PDEs (bijvoorbeeld de Navier Stokes vergelijkingen of gecombineerd warmte en vochttransport) kan simuleren. De gebruiker kan zich hierbij concentreren op het model (PDE coëfficiënten op het domein en op de rand) en hoeft maar weinig tijd te spenderen aan het verkrijgen en visualiseren van de oplossing. FemLab is een MatLab toolbox en bezit een export mogelijkheid van modellen naar Simulink. Deze combinatie kan heel interessant zijn voor de volgende mogelijke gebruikers: 1) Wetenschappers: Zij kunnen zich volledig concentreren op de fysica van de modellen en de validatie ervan. 2) Ingenieurs: Zij kunnen reeds ontwikkelde modellen gebruiken, aanpassen en combineren, voor het oplossen van typische ingenieursproblemen gebruikmakend van bijvoorbeeld Simulink. Doordat het relatief eenvoudig is om Grafische User Interfaces te bouwen voor modellen in Simulink, kunnen ingenieurs gebruiksvriendelijk applicaties ontwikkelen voor eindgebruikers. 3) Ontwerpers: Zij kunnen gebruik maken van de gebruiksvriendelijk applicaties waardoor ze zich alleen hoeven te concentreren op de invoergegevens en de simulatieresultaten.

2. FEMLAB FemLab [1] is een toolbox geschreven in MatLab [4]. Het simuleert systemen van gekoppelde PDEs (tot 32 onafhankelijke variabelen). De gespecificeerde PDEs mogen niet lineair zijn en tijdsafhankelijk en kunnen worden toegepast op een 1D, 2D of 3D geometrie. De geometrie, PDEs en randvoorwaarden zijn gedefinieerd door coëfficiënten in een bepaalde datastructuur. De PDEs en randvoorwaarden worden beschreven door een coëfficiënt vorm:

∂u − ∇ ⋅ ( c ⊗ ∇ u + α ⊗ u − γ ) + β ⊗ ∇ u + au = f ∂t n ⋅ ( c ⊗ ∇ u + α ⊗ u − γ ) + qu = g − h T λ , hu = r d

a

in Ω

(1)

on ∂ Ω

of een algemene vorm:

d

a

∂u + ∇ ⋅Γ = F ∂t

in Ω (2)

T

 ∂R  − n ⋅Γ = G +   λ  ∂u 

, 0 = R

on

∂Ω

met u is de onbekende oplossing, Ω is een begrensd domein, ∂Ω is de rand van het domein, n is de naar buiten gerichte normaalvector. Alle andere symbolen zijn coëfficiënten.

Een compleet voorbeeld Om te laten zien hoe FemLab werkt, wordt een simpel voorbeeld compleet uitgewerkt. In figuur 1 zijn de geometrie, PDE (zie (3)) met randvoorwaarden gegeven voor een eenvoudig koudebrug probleem:

b1 b2

concrete

b3

insulation b4

b7 b5

∇ ⋅ ( K∇ T ) = 0

b6

(3)

Boundary values, length b1: heat flux=he*(Te-T), 1m b2: heat flux=0, 0.2 m b3: heat flux=0, 0.1 m b4: heat flux=hi*(Ti-T), 0.8 m b5: heat flux=hi*(Ti-T), 0.7 m b6: heat flux=0, 0.2 m b7: heat flux=0, 0.8 m h : heat transfer coefficients T : Temperature i : internal e: external K: heat conduction coefficient

Figuur 1. De geometrie en model van een 2D koudebrug probleem. Gebruikmakend van de coëfficiënt vorm (1) en de PDE (3) met randvoorwaarden, betekent dit dat alle coëfficiënten van (1) nul zijn behalve: c={Kconcr , Kinsul} en q={hi, he }and g={hi*Ti, he*Te }.

In figuur 2 is de complete FemLab code gegeven voor het simuleren van bovenstaand probleem (default waarden van alle coëfficiënten van de PDE en randvoorwaarden zijn nul). Het automatisch gegenereerde grid is weergegeven in figuur 2 en de oplossing in figuur 3.

%CONSTANTS hi=7.7; %heat transfer coefficient internal he=25; %heat transfer coefficient external Ti=20; %air temperature internal Te=-10; %air temperature external Kconcr=1; %heat conduction concrete Kinsul=0.03; %heat conduction insulation %GEOMETRY: poly2(xdata,ydata) ; 2D polygon CONCR=poly2([0 0 1 1 0.2 0.2],[0 1 1 0.8 0.8 0]); %concrete INSUL=poly2([0.2 0.2 1 1],[0.7 0.8 0.8 0.7]); %insulation fem.geom=CONCR+INSUL; % fem geometry fem.dim=1; % One component (temperature) %COEFFICIENTS OF THE PDE/Boundary problem fem.equ.c={Kconcr Kinsul }; fem.bnd.g={he*Te 0 0 hi*Ti hi*Ti 0 0}; fem.bnd.q={he 0 0 hi hi 0 0}; fem.mesh=meshinit(fem); % intialize mesh fem.sol=femlin(fem); % solve linear, steady problem %OUTPUT MESH, SOLUTION meshplot(fem) q=posteval(fem,'u'); postplot(fem,'tridata',q,'tribar','on')

Figuur 2. De complete code (% betekent commentaar)

Figuur 3. Het grid.

Figuur 4. De oplossing (temperatuurverdeling in de koudebrug)

3. VALIDATIE VAN FEMLAB SIMULATIES MET METINGEN 3.1 1D vochttransport in een poreus materiaal Vochttransport in een poreus bouwmateriaal, beschreven door [5] kan gemodelleerd worden door:

∂θ ∂θ = ∇ ⋅ ( D (θ ) ⋅ ∇θ ) on 0 < x < 0.024 , θ = 0 .27 on x = 0, = 0 on x = 0 .024 ∂t ∂x

(4)

met θ is het vochtgehalte (m3.m-3) en D is de diffusiecoëfficiënt (m2s-1). In figuur 5 is de FemLab code gegeven voor het simuleren van bovenstaand model. De resultaten van gemeten waterabsorptie profielen voor verschillende materialen [5] is gegeven in figuur 6. Daarnaast zijn de simulatie resultaten weergegeven in figuur 7. De simulatieresultaten komen zeer goed overeen met de metingen.

fem.geom=solid1([0 0.024]); %geometry Dw=['(3.2e-9)*exp(29*u)']; %diffusivity fem.dim=1; %one component fem.equ.c={Dw}; %c equals diffusivity fem.equ.da={1}; %da equals 1 fem.bnd.h={1 0}; % boundary value fem.bnd.r={0.27 0}; % boundary value fem.mesh=meshinit(fem);fem.mesh=meshrefine(fem); fem.init={0}; fem.sol=femtime(fem,'tlist',[0:10:100]); Figuur 5. FemLab code voor het oplossen van een 1D diffusieprobleem

Figuur 6. Gemeten vochtprofielen versus lambda (= x.t-1/2)

Figuur 7. Gesimuleerde vochtprofielen versus lambda (= x.t-1/2)

3.2 2D Luchtstroming in een ruimte Dit voorbeeld, gebaseerd op [6], beschrijft de snelheid- , druk- en temperatuursverdeling in een kamer die verwarmd wordt door een stroom van hete lucht. Het probleem wordt gemodelleerd door 2D, incompressibele laminaire flow met Boussinesq benadering en met constante fysische eigenschappen. In figuur 8 zijn het model, de geometrie, de randvoorwaarden en het belangrijkste gedeelte van de code weergegeven. De coëfficiënten van de PDEs en randvoorwaarden zijn nu beschreven m.b.v. de algemene vorm (2). Continuty equation ( * are all omitted)

∂u ∂v + =0 ∂x ∂y

u* =

y u v x , v* = , x* = , y* = Uo Uo W1 W1

L* =

T − Tw L p , p* = , T* = 2 W1 ρUo T0 − Tw

U-momentum equation

∂ (uu ) ∂ (vu ) ∂p 1 2 + =− + ∇ u ∂x ∂y ∂x Re

V-momentum equation

∂ (uv) ∂ (vv) ∂p 1 2 Gr + =− + ∇ v+ T ∂x ∂y ∂y Re Re 2

Gr =

Energy equation

∂ (uT ) ∂ (vT ) 1 + = ∇ 2T ∂x ∂y Re Pr

Pr =

gβ∆TW1 v C pµ

2

3

, Re =

UoW1 , v

k

eta=1/Re;beta=Gr/(Re*Re);alpha=1/(Re*Pr); fem.equ.da={{1; 1; 0; 1}}; fem.equ.f={{'-(u1.*u1x+u2.*u1y+u3x)';... '-(u1.*u2x+u2.*u2y+u3y)+beta*u4';... '-(u1x+u2y)';... '-(u1.*u4x+u2.*u4y)'}}; fem.equ.ga={{{'-eta*u1x'; '-eta*u1y'};... {'-eta*u2x'; '-eta*u2y'};... 0;... {'-alpha*u4x'; '-alpha*u4y'}}}; fem=femdiff(fem); The boundary conditions are: At the left, right, top and bottom walls: u=0, v=0, T=0. At the inlet u=1, v=0, T=1. At the outlet Neuman conditions for u,v and T

Figuur 8. Model, geometrie, randvoorwaarden en FemLab code (u1=u*,u2=v*,u3=p*,u4=T*)

De resultaten van [6] en van de FemLab simulaties zijn vergeleken in de volgende figuur:

Figuur 9. Vergelijking van de gevalideerde simulatie resultaten van [6] (links) met de FemLab simulatie resultaten (rechts).

Uit figuur 9 blijkt dat ook voor dit sterk niet lineair proces de resultaten van FemLab goed overeen komen met literatuurwaarden.

4. GEAVANCEERDE APPLICATIES

4.1 2D Convectieve luchtstroming rondom een convector Het doel van dit project [8] is om de luchtstroming rondom een convector geplaatst onder een koud vlak te onderzoeken. De convector wordt toegepast als Lage Temperatuur Verwarmingssysteem (LTV) . De kernvraag is bij welke gemiddelde oppervlakte temperatuur van de convector de koudeval wordt gecompenseerd. Naast experimenten zijn overeenkomstige simulaties verricht met de volgende geometrie:

Figuur 10. De geometrie van de kamer (3m x 2.5m), de convector en het koud vlak Voor het verkrijgen van simulatieresultaten is het model uit de vorige paragraaf gebruikt met aanpassing van de randvoorwaarden. In figuur 11 en 12 zijn de na 2 minuten ontstane luchttemperaturen weergeven (in K), bij een convector oppervlaktetemperatuur van 25 oC en resp. 37.5 oC. De luchttemperatuur aan het begin bedraagt 20 oC.

Figuur 11. De luchttemperatuur (in K) na 2 minuten, bij een convector temperatuur van 25 oC

Figuur 12. De luchttemperatuur (in K) na 2 minuten bij een convectortemperatuur van 37.5 oC Figuur 11 en 12 laten zien dat een convector oppervlakte temperatuur van 25 oC te weinig is voor de koudeval compensatie in tegenstelling tot convector oppervlakte temperatuur van 37.5 oC, die wel voldoende is. Het model is gevalideerd met metingen [8] voor de gegeven geometrie en is gebruikt voor het simuleren van de luchtstroming bij verschillende temperaturen voor de koude wand, de convector en de begintemperatuur van de kamer.

4.2 Een FemLab model gekoppeld aan een gebouwmodel in Simulink Het doel van dit project is te onderzoek hoe een FemLab model (bijvoorbeeld het koudebrug model) kan worden gekoppeld aan een model in Simulink (bijvoorbeeld een gebouwmodel). Het model van figuur 4 is geëxporteerd naar Simulink, gebruikmakend van de standaard export optie van FemLab en gekoppeld in aan gebouwmodel [9]. In figuur 13 is het complete Simulink model weergegeven:

Figuur 13. De koppeling van een koudebrug model met een gebouwmodel in Simulink.

Het hierboven gepresenteerde model is gebruikt om na te gaan wat het effect van nachtverlaging is op de relatieve vochtigheid nabij het koudste gedeelte van de koudebrug in het vertrek. Hiertoe worden gedetailleerde temperaturen van de koudebrug gesimuleerd alsmede de extra warmtestroom van en naar de koudebrug als functie van de tijd. In figuur 14 zijn de resultaten weergeven voor de situatie zonder nachtverlaging en figuur 15 met nachtverlaging.

Figuur 14. Interne en externe luchttemperatuur, minimale en maximale oppervlakte temperatuur van de koudebrug en bijbehorende relatieve vochtigheid versus de tijd [dagen] met een constante zone temperatuur.

Figuur 15. Interne en externe luchttemperatuur, minimale en maximale oppervlakte temperatuur van de koudebrug en bijbehorende relatieve vochtigheid versus de tijd [dagen] met nachtverlaging. Bovenstaande grafieken laten duidelijk het effect zien van een nachtverlaging op de relatieve vochtigheid nabij een koudebrug. Verder heeft dit project aangetoond dat lineair modellen in FemLab gemakkelijk geëxporteerd kunnen worden naar Simulink. De simulatietijd blijft redelijk beperkt: Bovenstaande model heeft 3 minuten rekentijd nodig op een Pentium 600MHz.

4.3 3D gecombineerd warmte- en vochttransport Het doel van dit project is om te onderzoeken hoe 3D warmte- en vochttransport in poreuze materialen gesimuleerd kan worden in FemLab. Het probleem is gebaseerd op [7], maar is nu uitgebreid naar een 3D geometrie. Een proefstuk van steen, initieel droog en aan twee zijden gesealed, wordt op t=0 met de onderkant 1 cm onder water geplaatst. In figuur 16 zijn het model, de geometrie, de randvoorwaarden en de FemLab code (gebaseerd op de coëfficiënten vorm (1)), weergegeven. In figuur 17 zijn de materiaaleigenschappen weergegeven en in figuur 18 zijn de gesimuleerde 3D vochtprofielen afgebeeld.

∂w = ∇ ⋅ ( Dw ( w ) ∇ w + DF (T , w ) ∇ T ) ∂t ∂T = ∇ ⋅ ( DT ( w ) ∇ T ) ∂t w (t = 0) = w 0 , T (t = 0 ) = T 0 , sealed sides : ∇ T = ∇ w = 0 open sides : ∇ T = A (T − Ti ), ∇ w = B ( w − wi ) w ( z = 0 ) = w max , T ( z = 0 ) = T 0 Lx=0.11;Ly=0.07;Lz=0.49; Ti=23;wi=16;T0=20;w0=16; fem.geom=block3(Lx,Ly,Lz); fem.mesh=meshinit(fem); fem.dim={'w' 'T'}; fem.form='coefficient'; fem.equ.da={{1; 1}}; fem.equ.c={'Dw(w)' 'DF(w,T)' 0 DT(w)'}; fem.init={{wi T0}} fem.sol=femtime(fem);

Figuur 16. Model, geometrie, randvoorwaarden en FemLab code

Figuur 17. De materiaaleigenschappen voor de PDE coëfficiënten Dw,DT en DF

Figuur 18. Visualisatie van de vochtprofielen in 3D voor achtereenvolgens t=0, t=5 en t=70 dagen. Bovenstaande resultaten zijn indirect gevalideerd met metingen. Hiermee is het één van eerste modellen die in staat is om accuraat 3D warmte- en vochttransport te simuleren.

4.4 Een luchtstromingsmodel geregeld door een aan/uit regelaar in Simulink Het doel van dit project is het koppelen van een luchtstromingsmodel in FemLab met een regelaar in Simulink. Omdat het luchtstromingsmodel erg niet lineair is, kan de standaard export mogelijkheid van FemLab niet gebruikt worden. Dit probleem is opgelost door het schrijven van een discrete S-functie in Simulink. De S-functie simuleert een luchtstromingsprobleem met een bepaalde tijdstap (1 sec), met startwaarden uit de vorige tijdstap. Na iedere tijdstap wordt een gedeelte van de oplossing (lokatie van de sensor) geëxporteerd naar Simulink als input voor de aan/uit regelaar. Afhankelijk van de aan/uit regelaar, wordt al dan niet de randvoorwaarde voor de inblaaslucht aangepast. In het volgende voorbeeld is het model van figuur 9 gebruikt. De temperatuur van de inblaaslucht wordt nu geregeld door een aan/uit regelaar met hysterese (Relay) van Simulink. Als de temperatuur van de sensor (zie o in onderstaande figuur) boven de 20.5 oC komt, wordt er lucht ingeblazen van 17 oC en als de sensortemperatuur onder de 19.5 oC is dit 22 oC. De begintemperatuur van de ruimte is 17 oC. De volgende figuur laat het resultaat zien:

Figuur 20

Figuur 19. Het Simulink model, de temperatuur van de sensor en de output van de aan/uit regelaar versus de tijd (s),de temperatuurverdeling in de ruimte op t= 6 sec. (warme lucht wordt ingeblazen) en op t = 12 sec. (koude lucht wordt ingeblazen) (o = sensor positie).

De bovenstaande resultaten laten zien dat ook sterk niet lineaire modellen in FemLab, geëxporteerd kunnen worden naar Simulink, door het schrijven van een S-function.

5. TOEPASSINGSMOGELIJKHEDEN IN HET (VROEGE) ONTWERPSTADIUM De volgende aspecten zijn van belang voor het toepassen van FemLab/Simulink/Matlab in het (vroege) ontwerpstadium: • •

•

•

Betrouwbaarheid. Voor handhaving van de kwaliteit is het noodzakelijk een goede documentatie van modellen en validaties toe te voegen. Eenvoudig in gebruik. Complete modellen (of applicaties), waarvan men verwacht dat er in de praktijk veel behoefte aan is, kunnen worden voorzien van een gebruiksvriendelijke user-interface. Een 'beginnend' gebruiker kan hiermee zonder inhoudelijk kennis van de modellen, direct simulatieresultaten verkrijgen. Dit is een van de uitgangspunten van het Climasim project [10]. Flexibiliteit. Dit is bijvoorbeeld van belang als het ontwerp net even anders is dan de standaard oplossing die aanwezig is. Door de open structuur van de FemLab/Simulink/Matlab omgeving is een 'ervaren' gebruiker in staat om binnen een beperkte tijd, modellen aan te passen en/of eigen modellen toe te voegen. De nieuwe applicatie kan dan weer door een 'beginnend' gebruiker gebruikt worden. Een 'ervaren' gebruiker zou een ingenieur kunnen zijn die belast is met simulatievraagstukken binnen een bedrijf of overheidsinstelling. Compleetheid. M.b.v. deze simulatieomgeving, is een 'meer ervaren' gebruiker in staat de meest uiteenlopende problemen te modelleren op het gebied van bouwfysica, installatietechniek, regeltechniek en combinaties hiervan. Echter, het ontwikkelen van een geheel nieuw model vergt enige tijd, zodat deze niet meteen ingezet kan worden in het ontwerpproces.

6. CONCLUSIES FemLab is geëvalueerd voor een aantal problemen op het gebied van bouwfysica en installatietechniek die gebaseerd zijn op partiele differentiaal vergelijkingen (PDEs). Voorbeelden van typische bouwfysische problemen zoals dynamische luchtstroming en gecombineerd warmte- en vochttransport in poreuze materialen zijn m.b.v. FemLab gesimuleerd en gevalideerd. De simulatie resultaten hiervan zijn in goede overeenstemming met metingen. Om te laten zien hoe krachtig de FemLab software is, ook in combinatie met Simulink, worden de volgende applicaties behandeld: 1) het modeleren van 3D gecombineerd warmte- en vochttransport in een materiaal, 2) het modeleren van koudeval compensatie door een (put)convector, 3) het dynamisch simuleren van een 2D koudebrug in combinatie met een gebouwmodel, 4) het regelen van een luchtstromingsmodel in FemLab met een aan/uit regelaar van Simulink. FemLab en Simulink werken in de MatLab omgeving. Hierdoor is mogelijk om de reeds aanwezige visualisatie tools, toolboxen en alle andere programma's geschreven in Matlab te gebruiken. Mede dankzij de recente ontwikkelingen van Matlab, is de combinatie FemLab/Simulink/Matlab een krachtige en flexibele software omgeving voor het oplossen van wetenschappelijke en technische problemen op het gebied van bouwfysica en installatietechniek.

LITERATUUR [1] COMSOL AB, FEMLAB Version 2.0 pre, Reference Manual, September 2000 [2] PDE Solutions Inc, FlexPDE Manual, Version 2.14, August 1999 [3] Macsyma, Inc. PDEase2 version 2.5.4, Reference Manual, June 1995. [4] The Mathworks, Inc. MatLab Manual, Version 5.3. 1998 [5] Brocken, H., Moisture transport in brick masonry, Ph.D. thesis, Eindhoven University of Technology, December 1998 [6] Sinha, S.L., et al, Numerical simulation of two-dimensional room air flow with and without buoyancy, Energy and Buildings 32(2000) 121-129 [7] Kunzel, H.M., IEA Annex 24 HAMTIE, Final Report – Task 1 [8] Schijndel, A.W.M. van, Aarle, M.A.P. van, Investigation on the convection principle for convectors operating at low temperatures (Dutch), Report Eindhoven University of Technology, June 2001 [9] Verdonschot J.K.M. et al, Ontwikkeling van een simulatieomgeving voor klimaatsystemen, Proceedings IBPSA Symposium december 2001 Petten [10] Taal, A.C, Project Climasim, Stand van Zaken, Proceedings IBPSA Symposium december 2001 Petten