LORENE
hoteos.C
1 /*
2  * Methods of the class Hot_eos.
3  *
4  * (see file hoteos.h for documentation).
5  *
6  */
7 
8 /*
9  * Copyright (c) 2015 Jerome Novak
10  *
11  * This file is part of LORENE.
12  *
13  * LORENE is free software; you can redistribute it and/or modify
14  * it under the terms of the GNU General Public License version 2
15  * as published by the Free Software Foundation.
16  *
17  * LORENE is distributed in the hope that it will be useful,
18  * but WITHOUT ANY WARRANTY; without even the implied warranty of
19  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20  * GNU General Public License for more details.
21  *
22  * You should have received a copy of the GNU General Public License
23  * along with LORENE; if not, write to the Free Software
24  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
25  *
26  */
27 
28 char hoteos_C[] = "$Header: /cvsroot/Lorene/C++/Source/Eos/hoteos.C,v 1.2 2015/09/10 13:28:18 j_novak Exp $" ;
29 
30 /*
31  * $Id: hoteos.C,v 1.2 2015/09/10 13:28:18 j_novak Exp $
32  * $Log $
33  *
34  * $Header $
35  *
36  */
37 
38 // C headers
39 #include <cassert>
40 
41 // Lorene headers
42 #include "hoteos.h"
43 #include "eos.h"
44 #include "tensor.h"
45 #include "utilitaires.h"
46 
47 namespace Lorene {
48 
49  //--------------//
50  // Constructors //
51  //--------------//
52 
53 
54 // Standard constructor without name
55 // ---------------------------------
57 
58  // Standard constructor with name
59  // ---------------------------------
60  Hot_eos::Hot_eos(const string& name_i):name(name_i){
61  set_der_0x0() ;
62  }
63 
64  Hot_eos::Hot_eos(const char* name_i):name(name_i){
65  set_der_0x0() ;
66  }
67 
68  // Copy constructor
69  // ----------------
70  Hot_eos::Hot_eos(const Hot_eos& eos_i):name(eos_i.name){
71  set_der_0x0() ;
72  }
73 
74  // Constructor from a binary file
75  // ------------------------------
76  Hot_eos::Hot_eos(FILE* fich){
77 
78  int taille ;
79  fread(&taille, sizeof(int), 1, fich) ;
80  assert(taille > 0) ;
81  char* t_name = new char[taille] ;
82  fread(t_name, sizeof(char), taille, fich) ;
83  set_name(t_name) ;
84  delete [] t_name ;
85  set_der_0x0() ;
86  }
87 
88  // Constructor from a formatted file
89  // ---------------------------------
90  Hot_eos::Hot_eos(ifstream& fich){
91 
92  char t_name[100] ;
93  fich.getline(t_name, 100) ;
94  set_name(t_name) ;
95  set_der_0x0() ;
96  }
97 
98  //--------------//
99  // Destructor //
100  //--------------//
103  }
104 
105  //----------------------------------//
106  // Management of derived quantities //
107  //----------------------------------//
108 
109  void Hot_eos::del_deriv() const {
110  if (p_cold_eos != 0x0) delete p_cold_eos ;
111  set_der_0x0() ;
112  }
113 
114  void Hot_eos::set_der_0x0() const {
115  p_cold_eos = 0x0 ;
116  }
117 
118 void Hot_eos::set_name(const char* name_i) {
119 
120  name.assign(name_i) ;
121 
122 }
123 
124 
125  //------------//
126  // Outputs //
127  //------------//
128 
129 void Hot_eos::sauve(FILE* fich) const {
130 
131  int ident = identify() ;
132  fwrite_be(&ident, sizeof(int), 1, fich) ;
133 
134  int taille = int(name.size()) ;
135  fwrite_be(&taille, sizeof(int), 1, fich) ;
136  fwrite(name.c_str(), sizeof(char), name.size(), fich) ;
137 
138 }
139 
140 
141 
142 
143 ostream& operator<<(ostream& ost, const Hot_eos& eqetat) {
144  ost << eqetat.get_name() << endl ;
145  eqetat >> ost ;
146  return ost ;
147 }
148 
149  //-------------------------------//
150  // Generic computational routine //
151  //-------------------------------//
152 
153  void Hot_eos::calcule(const Scalar& ent, const Scalar& sb, int nzet, int l_min,
154  double (Hot_eos::*fait)(double, double) const,
155  Scalar& resu) const {
156 
157  assert(ent.get_etat() != ETATNONDEF) ;
158  assert(sb.get_etat() != ETATNONDEF) ;
159 
160  const Map* mp = &(ent.get_mp()) ; // Mapping
161 
162  const Mg3d* mg = mp->get_mg() ; // Multi-grid
163 
164  int nz = mg->get_nzone() ; // total number of domains
165 
166  if (ent.get_etat() == ETATZERO) {
167  resu.set_etat_zero() ;
168  return ;
169  }
170 
171  assert(ent.get_etat() == ETATQCQ) ;
172  const Valeur& vent = ent.get_spectral_va() ;
173  vent.coef_i() ; // the values in the configuration space are required
174 
175  const Valeur* vsb = &sb.get_spectral_va() ;
176  Valeur vzero(mg) ;
177  if (sb.get_etat() == ETATZERO) {
178  vzero.annule_hard() ;
179  vsb = &vzero ;
180  }
181 
182  assert(vsb->get_mg() == vent.get_mg()) ;
183 
184  // Preparations for a point by point computation:
185  resu.set_etat_qcq() ;
186  Valeur& vresu = resu.set_spectral_va() ;
187  vresu.set_etat_c_qcq() ;
188  vresu.c->set_etat_qcq() ;
189 
190  // Loop on domains where the computation has to be done :
191  for (int l = l_min; l< l_min + nzet; l++) {
192 
193  assert(l>=0) ;
194  assert(l<nz) ;
195 
196  bool tsb0 = false ;
197  Tbl* tent = vent.c->t[l] ;
198  Tbl* tsb = vsb->c->t[l] ;
199  Tbl* tresu = vresu.c->t[l] ;
200 
201  if (tent->get_etat() == ETATZERO) {
202  tresu->set_etat_zero() ;
203  }
204  else {
205  assert( tent->get_etat() == ETATQCQ ) ;
206  tresu->set_etat_qcq() ;
207  if (tsb->get_etat() == ETATZERO) {
208  tsb0 = true ;
209  tsb = new Tbl(tent->dim) ;
210  tsb->annule_hard() ;
211  }
212 
213  for (int i=0; i<tent->get_taille(); i++) {
214 
215  tresu->t[i] = (this->*fait)( tent->t[i], tsb->t[i] ) ;
216  }
217 
218  } // End of the case where ent != 0 in the considered domain
219  if (tsb0) delete tsb ;
220  } // End of the loop on domains where the computation had to be done
221 
222  // resu is set to zero in the other domains :
223 
224  if (l_min > 0) {
225  resu.annule(0, l_min-1) ;
226  }
227 
228  if (l_min + nzet < nz) {
229  resu.annule(l_min + nzet, nz - 1) ;
230  }
231  }
232 
233  // Baryon density from enthalpy
234  //------------------------------
235 
236  Scalar Hot_eos::nbar_Hs(const Scalar& ent, const Scalar& sb, int nzet, int l_min)
237  const {
238 
239  Scalar resu(ent.get_mp()) ;
240 
241  calcule(ent, sb, nzet, l_min, &Hot_eos::nbar_Hs_p, resu) ;
242 
243  return resu ;
244 
245  }
246 
247 
248 
249  // Energy density from enthalpy
250  //------------------------------
251 
252  Scalar Hot_eos::ener_Hs(const Scalar& ent, const Scalar& sb, int nzet, int l_min)
253  const {
254 
255  Scalar resu(ent.get_mp()) ;
256 
257  calcule(ent, sb, nzet, l_min, &Hot_eos::ener_Hs_p, resu) ;
258 
259  return resu ;
260 
261  }
262 
263  // Pressure from enthalpy
264  //-----------------------
265 
266  Scalar Hot_eos::press_Hs(const Scalar& ent, const Scalar& sb, int nzet, int l_min)
267  const {
268 
269  Scalar resu(ent.get_mp()) ;
270 
271  calcule(ent, sb, nzet, l_min, &Hot_eos::press_Hs_p, resu) ;
272 
273  return resu ;
274 
275  }
276 
277  // Temperature from enthalpy
278  //--------------------------
279 
280  Scalar Hot_eos::temp_Hs(const Scalar& ent, const Scalar& sb, int nzet, int l_min)
281  const {
282 
283  Scalar resu(ent.get_mp()) ;
284 
285  calcule(ent, sb, nzet, l_min, &Hot_eos::temp_Hs_p, resu) ;
286 
287  return resu ;
288 
289  }
290 
291 }
Lorene::Hot_eos::nbar_Hs_p
virtual double nbar_Hs_p(double ent, double sb) const =0
Computes the baryon density from the log-enthalpy and entropy per baryon (virtual function implemente...
Lorene::Valeur
Values and coefficients of a (real-value) function.
Definition: valeur.h:287
Lorene::Valeur::set_etat_c_qcq
void set_etat_c_qcq()
Sets the logical state to ETATQCQ (ordinary state) for values in the configuration space (Mtbl c ).
Definition: valeur.C:701
Lorene::Mg3d
Multi-domain grid.
Definition: grilles.h:273
Lorene::Hot_eos::temp_Hs
Scalar temp_Hs(const Scalar &ent, const Scalar &sb, int nzet, int l_min=0) const
Computes the temperature field from the log-enthalpy field and entropy per baryon.
Definition: hoteos.C:280
Lorene
Lorene prototypes.
Definition: app_hor.h:64
Lorene::Tbl::annule_hard
void annule_hard()
Sets the Tbl to zero in a hard way.
Definition: tbl.C:372
Lorene::Map::get_mg
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:765
Lorene::Hot_eos::Hot_eos
Hot_eos()
Standard constructor.
Definition: hoteos.C:56
Lorene::Valeur::get_mg
const Mg3d * get_mg() const
Returns the Mg3d on which the this is defined.
Definition: valeur.h:729
Lorene::Scalar
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:387
Lorene::Tensor::get_mp
const Map & get_mp() const
Returns the mapping.
Definition: tensor.h:861
Lorene::Tbl::dim
Dim_tbl dim
Number of dimensions, size,...
Definition: tbl.h:172
Lorene::Valeur::c
Mtbl * c
Values of the function at the points of the multi-grid
Definition: valeur.h:299
Lorene::Tbl::t
double * t
The array of double.
Definition: tbl.h:173
Lorene::Scalar::get_spectral_va
const Valeur & get_spectral_va() const
Returns va (read only version)
Definition: scalar.h:601
Lorene::Scalar::set_spectral_va
Valeur & set_spectral_va()
Returns va (read/write version)
Definition: scalar.h:604
Lorene::Tbl::set_etat_qcq
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tbl.C:361
Lorene::Hot_eos::ener_Hs
Scalar ener_Hs(const Scalar &ent, const Scalar &sb, int nzet, int l_min=0) const
Computes the total energy density from the log-enthalpy and entropy per baryon.
Definition: hoteos.C:252
Lorene::Valeur::annule_hard
void annule_hard()
Sets the Valeur to zero in a hard way.
Definition: valeur.C:723
Lorene::Hot_eos::nbar_Hs
Scalar nbar_Hs(const Scalar &ent, const Scalar &sb, int nzet, int l_min=0) const
Computes the baryon density field from the log-enthalpy field and entropy per baryon.
Definition: hoteos.C:236
Lorene::Hot_eos::identify
virtual int identify() const =0
Returns a number to identify the sub-classe of Hot_eos the object belongs to.
Lorene::Tbl::set_etat_zero
void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: tbl.C:347
Lorene::Tbl
Basic array class.
Definition: tbl.h:161
Lorene::Mg3d::get_nzone
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:448
Lorene::Scalar::set_etat_qcq
virtual void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: scalar.C:353
Lorene::Scalar::annule
virtual void annule(int l_min, int l_max)
Sets the Scalar to zero in several domains.
Definition: scalar.C:391
Lorene::Scalar::set_etat_zero
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: scalar.C:324
Lorene::Hot_eos::set_name
void set_name(const char *)
Sets the hot EOS name.
Definition: hoteos.C:118
Lorene::Hot_eos::ener_Hs_p
virtual double ener_Hs_p(double ent, double sb) const =0
Computes the total energy density from the log-enthalpy and entropy per baryon (virtual function impl...
Lorene::Hot_eos::press_Hs
Scalar press_Hs(const Scalar &ent, const Scalar &sb, int nzet, int l_min=0) const
Computes the pressure from the log-enthalpy and entropy per baryon.
Definition: hoteos.C:266
Lorene::Scalar::get_etat
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
Definition: scalar.h:554
Lorene::Hot_eos::press_Hs_p
virtual double press_Hs_p(double ent, double sb) const =0
Computes the pressure from the log-enthalpy and entropy per baryon (virtual function implemented in t...
Lorene::Hot_eos::~Hot_eos
virtual ~Hot_eos()
Destructor.
Definition: hoteos.C:101
Lorene::Hot_eos::sauve
virtual void sauve(FILE *) const
Save in a file.
Definition: hoteos.C:129
Lorene::Mtbl::t
Tbl ** t
Array (size nzone ) of pointers on the Tbl 's.
Definition: mtbl.h:132
Lorene::Hot_eos::del_deriv
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: hoteos.C:109
Lorene::fwrite_be
int fwrite_be(const int *aa, int size, int nb, FILE *fich)
Writes integer(s) into a binary file according to the big endian convention.
Definition: fwrite_be.C:70
Lorene::Hot_eos::calcule
void calcule(const Scalar &thermo1, const Scalar &thermo2, int nzet, int l_min, double(Hot_eos::*fait)(double, double) const, Scalar &resu) const
General computational method for Scalar 's.
Definition: hoteos.C:153
Lorene::Valeur::coef_i
void coef_i() const
Computes the physical value of *this.
Definition: valeur_coef_i.C:140
Lorene::Hot_eos::p_cold_eos
Eos * p_cold_eos
Corresponding cold Eos.
Definition: hoteos.h:108
Lorene::Tbl::get_taille
int get_taille() const
Gives the total size (ie dim.taille)
Definition: tbl.h:397
Lorene::Map
Base class for coordinate mappings.
Definition: map.h:670
Lorene::Hot_eos::get_name
const string & get_name() const
Returns the hot EOS name.
Definition: hoteos.h:120
Lorene::Hot_eos::name
string name
EOS name.
Definition: hoteos.h:72
Lorene::Tbl::get_etat
int get_etat() const
Gives the logical state.
Definition: tbl.h:394
Lorene::Hot_eos
Base class for temperature-dependent equations of state (abstract class).
Definition: hoteos.h:67
Lorene::Hot_eos::temp_Hs_p
virtual double temp_Hs_p(double ent, double sb) const =0
Computes the temperature from the log-enthalpy and entropy per baryon (virtual function implemented i...
Lorene::Mtbl::set_etat_qcq
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: mtbl.C:299
Lorene::Hot_eos::set_der_0x0
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: hoteos.C:114