Benutzer:Widakora/Grobmodell

aus Wikipedia, der freien Enzyklopädie
Zur Navigation springen Zur Suche springen
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 fr 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 fr 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 fr den Turmwirkungsgrad, s.u.}
{sekundaere Wetter-Werte}


{A:=ATM(TA,RF,Glob,GRECH); wegen Errechnen des Sonnenstandes in Grech zu aufwendig
fr 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 eingefhrt mit der Absicht,
die bekannte Gasgleichung fr Luft auch fr 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 fr 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 fr 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 Ausdrcken 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„rmebergangskoeffizient 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
natrlich mit dem im Kollektor bereinstimmen muss !
Gesucht ist also dasjenige dT, das fr š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 | Quelltext bearbeiten]
  1. 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.