Program Grobmodell; {widakora 8.12.2010}
(*dieses Programm ist unter TURBO PASCAL 5.0 lauffähig.
Es erhebt keinerlei Anspruch auf Programmier-Ästhetik,
Fehlerfreiheit oder Aktualität. Ich habe seit 10 Jahren
nicht mehr programmiert, zu Zeiten dieser TURBO-Version
war ich blutiger EDV- Anfänger mit trial-and-error-Programmierstil
Die Kosten-Angaben im Programm sind noch in DM und dem BMFT- Forschungsbericht[1] entnommen. Neue EUR- Angaben finden sich im Hauptartikel fürs Aufwindkraftwerk, der sich derzeit noch auf der Benutzer:widakora/Baustelle befindet
Dieses Grobmodell ist gegenüber dem 3-Faktoren-Modell insoweit eine Verbesserung, als hier die Kollektorverluste als Funktion vom Temperaturhub dT näherungsweise durch empirische Ansätze berechnet werden können, so dass die Aufwindgeschwindigkeit berechnet werden kann und im Einlaufbereich natürlich aus Kontinuitäts- und Massenerhaltungsgründen mit der Kollektor- Endgeschwindigkeit in Einklang gebracht werden muss, was natürlich wiederum Rückwirkungen auf den Kollektorwirkungsgrad und das Aufwind-erzeugende dT hat. Immerhin kann- immer im Zusammenspiel mit den korrespondierenden Baukosten- mit dem Turmradius als zusätzlichem Abmessungsparamter "gespielt" und eine physikalisch-bautechnische Groboptimierung durchgeführt werden. Was dieses Primitiv- Modell nicht leistet: eine dynamische Brücksichtigung der Wärmeleitung in (und aus) den(m) Boden, eine Berechnung der Reibungsdruckverluste (diese sind durch einen kumulativen empirischen Anlagen-cw-Wert berücksichtigt unter der Vorgabe, dass der Aufwind aus Erfahrung 15m/s nicht wesentlich übersteigen sollte, da sonst die Kollektorreibung zu groß oder die Kollektorkosten zu hoch werden wegen dann sehr großer Bauhöhen des Treibhauses. Ferner ist hier natürlich nur eine statische Momentaufnahme möglich z.B. für Peakleistung oder mittlere Jahresleistungen- in einem dynamischen Modell müssen alle Wetterparameter in Zeitintervallen von maximal 10min vorgegeben werden. Ausfürliche Modelle für die verschiedenen Anlagenbereiche finden sich z.B. in den Proceedings der "SCPT 2010 Solar Chimney Power Technology International Conference, Ruhr-Universität Bochum, Bergische Universität Wuppertal 28.-30. September 2010 *)
var Zaehler: integer; Pel,eta_turm,eta_C, S,Tvakelv,takelv, Tabd : extended; dT,vT,sa,rhowa:extended; R,E,K,Qbod,Vd : extended; rlw,alphaAbd,rhoAbd:extended; HSTort,A,Ta,Glob,relFeu,Pa,Va,wdpa:extended; Turm_Flaeche,Endradius_Vordach,Koll_Flaeche, AWK_Flaeche:extended; tvikelv, rhoi, poben , Dp,Dpturb :extended; QLuft,mpunkt_Turm,mpunkt_Koll :extended; Baukosten,Turmkosten,Maschinenkosten,Vordachkosten:extended; spezifische_Investitionskosten:extended; c0,help,help1,help2,help3,help4:extended; fertig:boolean;
{Naturkonstanten}
const
T0 = 273.15;
Sigma = 5.67E-8; { Stefan-Boltzmannkonstante (W/K^4/m^2) }
CpLuft = 1005.0; { Spez. Waermekap. der Luft (J/KG/K)=Ws/kg/K }
CpWas = 1848.0; { wasdampf }
C = 0.0341752; { Konstante C = g/Rl (K/M) }
Gama = 0.0097; { Adiabatenkoeffizient (K/M) }
RW = 461.5; { Gaskonst Wasser}
Rl = 287.05;
RldivRw = 0.6219935;
g = 9.81; {Schwerebeschleunigung m/s^2 }
pi = 3.14159265;
RAD = 0.017453293; {PI/180}
{Anlagen- Geometrie und Materialwerte}
EPSLW = 0.95; {lw Absorption Abdeckung}
TAULW = 0.001;{lw Transmission Abdeckung }
ABSKW = 0.06; {kw Absorption Abdeckung }
TAUKW = 0.90; {kw Transmission Abdeckung ,-1% bis -2% fuer Unterkonstr. }
EPSLWBOD = 0.93; {lw Emissionsgrad Boden }
ABSKWBOD = 0.90; {kw Absorptionsgrad " }
CPAB = 481.0;{Spezif. Waermekapazitaet Abdeckung (J/kg/K)}
RHOAB = 2580.0;{Dichte Abdeckung (kg/qbm)}
RadiusTurm = 96; {m}
RadiusVordach = 2400;{m}
HTurm = 1000;{m}
Xt = 0.84;
eta_masch = 0.721;
cw = 1.1;
spMk = 1.850; {spezifische Maschinen-Kosten in DM/Wpk}
spTK = 400; {spezifische Turm-Kosten in DM/qm Turmoberfl„che}
spVDK = 38; {spezifische Vordach-Kosten in DM/qm Kollektorfl„che}
{empirische Anlagen-Parameter aus Manzanares}
kappa = 0.47; {Koeffizient fr die N„herung der aktuellen Abdeckungstemperatur als
Funktion der Lufttemperatur und des mittleren Temperaturhubes dT nach
Gl. 2.3 S.96 nach []}
beta = 2.1; {Koeffizient fr die N„herung der aktuellen Abdeckungstemperatur als
Funktion der Boden-Oberfl„chentemperatur und des mittleren Tempera-turhubes dT nach
Gl. 2.4 S.96 nach []}
FUNCTION XHOCHY ( X,Y: EXTENDED) : EXTENDED;
VAR CH : CHAR;
{ Funktion potenziert die Basis x mit dem Exponenten y. }
begin
if x = 0 then
if y = 0 then
xhochy:=1
else
if y > 0 then xhochy:=0 else write('Fehler in xhochy ')
else
if x > 0 then
xhochy := exp ( y * ln (x) )
else
if y = trunc(y) then
if odd(trunc(y)) then
xhochy:=-xhochy(-x,y)
else
xhochy:= xhochy(-x,y)
else
BEGIN
writeln;
WRITELN('Fehler beim Aufruf von XhochY.....');
WRITELN('X-WERT : ',X:8:2,'Y- WERT : ',Y:8:2);
readln;
END; end;
FUNCTION WASSERDAMPFDRUCK(RF,TA:DOUBLE):DOUBLE;
{ Wasserdampfpartialdruck in Pascal }
{ Rel.Feuchte in %, Temp. in Grad C }
CONST RF1=6.107;
RF2=7.433;
RF3=234.9;
BEGIN WASSERDAMPFDRUCK:=RF1*RF*XHOCHY(10,RF2*(TA)/(RF3+TA)); END;
Function POLY( GAM, H, TKELV : extended ) : extended;
{ HOEHENFORMEL FUER POLYTROPE SCHICHTUNG P(H0-H) = P(H0)*F(GAM,H,T) }
begin
Poly := XhochY (1.0 + GAM*H/TKELV,C/GAM );
end;
begin{main}
fertig:=false;
{Standort-Wetter}
Glob :=1005;{1050;}{292;}{652;} {W/qm, Max-arithm. Mittel-LeistungsWichtg
der j„hrlichen Globalstrahlung}
Pa :=100000; {Pascal}
Va :=2.3; {m/s, Leistungsgewichtetes Mittel der Windgeschwindigkeit in 10m}
relFeu:=45; {relative Feuchte in %}
Ta :=21.1; {C, Leistungsgewichtetes Mittel der Aussentemperatur am Boden 2m in Celsius}
A := 371.5; {W/qm Leistungsgewichtetes Mittel der atmosph„rischen Gegenstrahlung}
HSTort:=400; {H”he Standort BARSTOW}
S :=1.017; {Korrekturfaktor fr den Turmwirkungsgrad, s.u.}
{sekundaere Wetter-Werte}
{A:=ATM(TA,RF,Glob,GRECH); wegen Errechnen des Sonnenstandes in Grech zu aufwendig
fr das Grobmodell, es w„ren eine Anzahl weiterer Spezifikationen wie
Klimatyp, Niederschlag, Geogr. Breite und H”he etc. notwendig}
Takelv := ta+t0; {K}
(* Reduzierung Pa auf geo.hoehe:50.3 186.6*)
PA := PA*xhochy((1-gama*HSTort/TAKELV),C/gama);
wdpa := Wasserdampfdruck(RelFeu,ta); {Wasserdampfpartialdruck in Pa}
rhowa:= wdpa/takelv/rw; {Dichte Wasserdampf}
sa := rldivrw*wdpa/(pa + wdpa*(rldivrw-1));
tvakelv := takelv*(1+(rw/rl-1)*sa);
{Die virtuelle Temperatur Tvakelv wird eingefhrt mit der Absicht,
die bekannte Gasgleichung fr Luft auch fr Wasserdampfgemische anwendbar zu machen}
{sekundaere Anlagen-Werte}
Endradius_Vordach:=(* 2 *) RadiusTurm;
Turm_Flaeche := Pi*xhochy(RadiusTurm,2);
AWK_Flaeche := Pi *xhochy(Radiusvordach,2);
Koll_Flaeche := AWK_Flaeche-turm_flaeche;
rlw := 1 - epslw - taulw; {Reflexion langwellig als grauer Strahler}
Zaehler:=1;
dT:=20; {Anfangswert fr Schritt- Iteration}
eta_turm:=S*g*HTurm /(Cpluft*Tvakelv); {Turmwirkungsgrad Gl. 2.53 S.133}
{ Juni Juli August Sept}
{S 1.076 1.057 1.031 1.017 }
{Leistungsgewichteter gemessener Monatsmittelwert des
Korrekturfaktors S fr die adiabat angen„hert errechnete Gesamtdruckdifferenz}
R:=Glob*(1-taukw-abskw + (1-abskwbod)*taukw*taukw);
{Bilanz kurzwellig (Reflexion) nach Gl. 2.1 S.92 []}
{Schleife mit dT- abh„ngigen Ausdrcken zum Massendurchsatz-Angleich}
repeat Tabd:=Takelv+kappa*dT;
{mittlere Folienoberfl„chentemperatur als Funktion des Temperaturhubs im Kollektor
und der Umgebungstemperatur in 2m H”he}
E:=sigma*(epslw*xhochy(Tabd,4)+(1-epslw-rlw)*epslwbod*xhochy(Ta+beta*dT,4))-A;
{Bilanz langwellig (IR- Emission) nach Gl. 2.2 S.96 []}
alphaAbd:= 5.4+3*Va; {Gl. 2.9 S.98 []}
{W„rmebergangskoeffizient Abdeckungsoberfl„che-Umgebung in W/qm/K}
K:=alphaAbd*kappa*dT;
{Konvektionverluste Anlage-Umgebung nach Gl. 2.13 S. 100 []}
VD :=0;
Qbod:=xhochy(Glob,2.2)*1.98e-5; (*empirisch: bei Boden wie in Manzanares geht während
Sonnen- und deshalb Leistungsintensiven Betriebszeiten
etwa 1/3 der Globalstrahlungsleistung über in die
Boden- Wärmespeicherung - dieser Anteil sorgt dann für die
Nachtleistung*)
(*qbod:=0; für Langzeit-Mittel ohne nennenswerten Nachtbetrieb*)
eta_C:=1 - (R+E+K+VD+Qbod)/Glob; {Gl. 2.43 S.125 []}
{Der Kollekorwirkungsgrad beschreibt, welcher Anteil der Globalstrahlung
auf die gesamte Kollektorfl„che in die W„rmeleistung der QLuft übergeht:}
QLuft:= eta_C*Koll_Flaeche*Glob; {Watt}
(*Diese Wärmeleistung führt je nach Luftgeschwindigkeit im Kollektor zu
einem gewissen Temperaturhub nach Qluft:=mpunkt*cpluft*dT;{Watt=kg/s*J=Ws/kg/K*K}
Gleichzeitig erzeugt der Temperaturhub am Turmfuá eine Druckdifferenz, die
je nach Xt-Einstellung der Turbinenbl„tter (oder auch Drosselklappen) in dynamischen Druck
gewandelt wird und so einen ganz bestimmten Massendurchsatz durch den Turm bedingt, der
natrlich mit dem im Kollektor bereinstimmen muss !
Gesucht ist also dasjenige dT, das fr šbereinstimmung der Massendurchs„tze sorgt,
oder bildlich gesprochen: dass die Luftgeschwindigkeit am Kollektorende mit der
am Turmfuá bereinstimmt*)
mpunkt_Koll:= QLuft/cpluft/dT; {Massenfluss im Kollektor}
tvikelv := (takelv+dT)*(1+(rw/rl-1)*sa);
rhoi := pa/Rl/tvikelv;
poben := pa*xhochy((1-gama* HTurm /tvakelv),c/gama);
Dp := pa - poben*poly(gama,HTurm ,tvikelv-gama*HTurm ) ;
Dpturb := Xt * Dp;
VT := SQRT(2*(Dp-Dpturb)/rhoi)/SQRT(1+CW);
{Massenstrom im Turm:}
mpunkt_turm:= rhoi*Turm_Flaeche*VT;
if mpunkt_turm>mpunkt_koll then dT:=dT*0.995 else dT:=dT*1.005; (* DIES ist der Punkt, an dem spielfreudige alternde Physiker und spielfreudige Programmier-Sauriers eingeladen sind, statt der uneleganten Schrittiteration z.B. die Newton-Rgel oder die Regula Falsi zu programmieren- reduziert die Rechnedurchgänge evtl. von ca. 350 auf nur 3-4 *) inc(Zaehler);
fertig:=abs(1- mpunkt_turm/mpunkt_koll)<0.0001;
until fertig or (zaehler>10000); (****************************einige Formeln zur Erläuterung: *****************************
Pel:=eta_Masch*dPTurb*vt*Turm_flaeche;{Gl.2.66 S.164 [3]}
Pel:=eta_masch*Xt*S*g*HTurm /Tvakelv*rhoi*Vt*dT*Turm_Flaeche;{Gl.2.67 S.164 [3]}
Pel:=Xt*eta_masch*eta_c*eta_turm*Glob*Koll_Flaeche; {Gl. 2.82 S. 178 [3]}
Turmkosten :=spTK*2*PI*RadiusTurm*HTurm ; {DM-- neue EUR- Kosten im Artikel über Aufwindkraftwerke !}
Maschinenkosten:=spMK*Pel; {DM}
Vordachkosten :=spVDK*Koll_Flaeche; {DM}
Baukosten:=Turmkosten+Maschinenkosten+Vordachkosten;
spezifische_Investitionskosten:=(Maschinenkosten+Turmkosten+Vordachkosten)/Pel; {DM/Wpk}
spezifische_Investitionskosten:= spMK+Turmkosten/Pel+Vordachkosten/Pel;
Pel:=Xt*eta_masch*eta_c*eta_turm*Glob*Koll_Flaeche; {Wpk}
{einsetzen Turmwirkungsgrad eta_turm= S*g*HTurm /(Cpluft*Tvakelv) dimensionslos}
Pel:=Xt*eta_masch*eta_c*S*g*HTurm /(Cpluft*Tvakelv)*Glob*Koll_Flaeche; {Wpk}
Ende einige Formeln zur Erläuterung
************************************************************************************)
{ersetzen eta_turm}
c0:=Xt*eta_masch*eta_c*S*g/(Cpluft*Tvakelv)*Glob;
Pel:=c0*HTurm *Koll_Flaeche;
(*mit
Xt = 0.84 {Anteil der Gesamt-Druckentnahme durch die Turbine }
eta_masch = 0.721 {Maschinenwirkungsgrad Getriebe,Generator etc. }
eta_c = 0.536 {KoLlektorwirkungsgrad }
S = 1.017 {Berücksichtigung der Abweich von adiabat Aussenschichtung}
CpLuft = 1005.0; {Spez. Waermekap. der Luft (J=Ws/KG/K) }
g = 9.81; {Schwerebeschleunigung m/s^2 }
Tvakelv = 295.57 {Aussentemperatur in 2m H”he (K) }
Glob = 1050 {Globalstrahlung W/qm }
ist c0 = 0.0033318 {Einheit 1/m }*)
{kürzen: Pel 1. Glied Hturm 2. Glied und Koll_Flaeche 3. Glied}
spezifische_Investitionskosten:=spMk+spTK*2*PI*RadiusTurm/(c0*Koll_Flaeche)+spVDK/(c0*Hturm);
(*198 MWpk MaschKosten= 366.3Mio DM Turmkosten=241.3 Mio Vordachkosten=686.53 Mio bei 198 Mpk*) (*=6.54 TDM/kWpk*) end.
Ein Vergleich der Ergebnisse mit denjenigen aus dem ausführlichen Modell (s.u. Anhang V, Tabel-le für die BARSTOW- Anlage) liefert:
Ausführl. Modell Grobmodell
Jahres-Peakleistung bei G=1050 W/qm: 194.8 MW 198 MW Mittlere Leistung bei G=292 W/qm: 51.7 MW 44.2 MW
Die Unterschiede rühren vor allem vom Ansatz für Qbod - im Grobmodell kann die tatsächliche Dynamik der Anlage mit Bodenspeicherung und -wiederabgabe natürlich nicht wiedergegeben werden.
Einzelnachweise
Bearbeiten- ↑ Schlaich J., Bergermann R., Friedrich, K., Haaf, W., Lautenschlager, H. 1986. Baureife Planung und Bau einer Demonstrationsanlage eines atmosphärenthermischen Aufwindkraftwerkes. Anwendungsnahe Auslegung größerer Einheiten und erweitertes Meßprogramm BMFT-FB-T 86-208 Stuttgart. Referenzfehler: Ungültiges
<ref>-Tag. Der Name „BMFT-FB-T 86-208“ wurde mehrere Male mit einem unterschiedlichen Inhalt definiert.