LORENE
eos.C
1/*
2 * Methods of class Eos.
3 *
4 * (see file eos.h for documentation).
5 */
6
7/*
8 * Copyright (c) 2000-2001 Eric Gourgoulhon
9 *
10 * This file is part of LORENE.
11 *
12 * LORENE is free software; you can redistribute it and/or modify
13 * it under the terms of the GNU General Public License as published by
14 * the Free Software Foundation; either version 2 of the License, or
15 * (at your option) any later version.
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
29char eos_C[] = "$Header: /cvsroot/Lorene/C++/Source/Eos/eos.C,v 1.7 2014/10/13 08:52:51 j_novak Exp $" ;
30
31/*
32 * $Id: eos.C,v 1.7 2014/10/13 08:52:51 j_novak Exp $
33 * $Log: eos.C,v $
34 * Revision 1.7 2014/10/13 08:52:51 j_novak
35 * Lorene classes and functions now belong to the namespace Lorene.
36 *
37 * Revision 1.6 2014/10/06 15:13:06 j_novak
38 * Modified #include directives to use c++ syntax.
39 *
40 * Revision 1.5 2004/01/14 15:59:42 f_limousin
41 * Add methos calcule, nbar_ent, der_bar_ent, der_press_ent and press_ent
42 * for Scalar's.
43 *
44 * Revision 1.4 2002/10/16 14:36:34 j_novak
45 * Reorganization of #include instructions of standard C++, in order to
46 * use experimental version 3 of gcc.
47 *
48 * Revision 1.3 2002/04/09 14:32:15 e_gourgoulhon
49 * 1/ Added extra parameters in EOS computational functions (argument par)
50 * 2/ New class MEos for multi-domain EOS
51 *
52 * Revision 1.2 2001/12/04 21:27:53 e_gourgoulhon
53 *
54 * All writing/reading to a binary file are now performed according to
55 * the big endian convention, whatever the system is big endian or
56 * small endian, thanks to the functions fwrite_be and fread_be
57 *
58 * Revision 1.1.1.1 2001/11/20 15:19:27 e_gourgoulhon
59 * LORENE
60 *
61 * Revision 2.5 2001/10/10 13:45:06 eric
62 * Modif Joachim : &(Eos::der_press_ent_p) -> &Eos::der_press_ent_p, etc...
63 * pour conformite au compilateur HP.
64 *
65 * Revision 2.4 2001/02/07 09:50:49 eric
66 * Suppression de la fonction derent_ent_p.
67 * Ajout des fonctions donnant les derivees de l'EOS:
68 * der_nbar_ent_p
69 * der_ener_ent_p
70 * der_press_ent_p
71 * ainsi que des fonctions Cmp associees.
72 *
73 * Revision 2.3 2000/02/14 14:33:05 eric
74 * Ajout des constructeurs par lecture de fichier formate.
75 *
76 * Revision 2.2 2000/01/21 15:17:28 eric
77 * fonction sauve: on ecrit en premier l'identificateur.
78 *
79 * Revision 2.1 2000/01/18 13:47:07 eric
80 * Premiere version operationnelle
81 *
82 * Revision 2.0 2000/01/18 10:46:15 eric
83 * /
84 *
85 *
86 * $Header: /cvsroot/Lorene/C++/Source/Eos/eos.C,v 1.7 2014/10/13 08:52:51 j_novak Exp $
87 *
88 */
89
90// Headers C
91#include <cstdlib>
92#include <cstring>
93
94// Headers Lorene
95#include "eos.h"
96#include "cmp.h"
97#include "scalar.h"
98#include "utilitaires.h"
99
100 //--------------//
101 // Constructors //
102 //--------------//
103
104
105// Standard constructor without name
106// ---------------------------------
107namespace Lorene {
109
110 set_name("") ;
111
112}
113
114// Standard constructor with name
115// ---------------------------------
116Eos::Eos(const char* name_i){
117
118 set_name(name_i) ;
119
120}
121
122// Copy constructor
123// ----------------
124Eos::Eos(const Eos& eos_i){
125
126 set_name(eos_i.name) ;
127
128}
129
130// Constructor from a binary file
131// ------------------------------
132Eos::Eos(FILE* fich){
133
134 fread(name, sizeof(char), 100, fich) ;
135
136}
137
138// Constructor from a formatted file
139// ---------------------------------
140Eos::Eos(ifstream& fich){
141
142 fich.getline(name, 100) ;
143
144}
145
146
147
148 //--------------//
149 // Destructor //
150 //--------------//
151
153
154 // does nothing
155
156}
157
158 //-------------------------//
159 // Manipulation of name //
160 //-------------------------//
161
162
163void Eos::set_name(const char* name_i) {
164
165 strncpy(name, name_i, 100) ;
166
167}
168
169const char* Eos::get_name() const {
170
171 return name ;
172
173}
174
175 //------------//
176 // Outputs //
177 //------------//
178
179void Eos::sauve(FILE* fich) const {
180
181 int ident = identify() ;
182 fwrite_be(&ident, sizeof(int), 1, fich) ;
183
184 fwrite(name, sizeof(char), 100, fich) ;
185
186}
187
188
189
190
191ostream& operator<<(ostream& ost, const Eos& eqetat) {
192 ost << eqetat.get_name() << endl ;
193 eqetat >> ost ;
194 return ost ;
195}
196
197
198 //-------------------------------//
199 // Generic computational routine //
200 //-------------------------------//
201
202
203void Eos::calcule(const Cmp& ent, int nzet, int l_min,
204 double (Eos::*fait)(double, const Param*) const, const Param* par, Cmp& resu) const {
205
206 assert(ent.get_etat() != ETATNONDEF) ;
207
208 const Map* mp = ent.get_mp() ; // Mapping
209
210
211 if (ent.get_etat() == ETATZERO) {
212 resu.set_etat_zero() ;
213 return ;
214 }
215
216 assert(ent.get_etat() == ETATQCQ) ;
217 const Valeur& vent = ent.va ;
218 vent.coef_i() ; // the values in the configuration space are required
219
220 const Mg3d* mg = mp->get_mg() ; // Multi-grid
221
222 int nz = mg->get_nzone() ; // total number of domains
223
224 // Preparations for a point by point computation:
225 resu.set_etat_qcq() ;
226 Valeur& vresu = resu.va ;
227 vresu.set_etat_c_qcq() ;
228 vresu.c->set_etat_qcq() ;
229
230 // Loop on domains where the computation has to be done :
231 for (int l = l_min; l< l_min + nzet; l++) {
232
233 assert(l>=0) ;
234 assert(l<nz) ;
235
236 Tbl* tent = vent.c->t[l] ;
237 Tbl* tresu = vresu.c->t[l] ;
238
239 if (tent->get_etat() == ETATZERO) {
240 tresu->set_etat_zero() ;
241 }
242 else {
243 assert( tent->get_etat() == ETATQCQ ) ;
244 tresu->set_etat_qcq() ;
245
246 for (int i=0; i<tent->get_taille(); i++) {
247
248 tresu->t[i] = (this->*fait)( tent->t[i], par ) ;
249 }
250
251 } // End of the case where ent != 0 in the considered domain
252
253 } // End of the loop on domains where the computation had to be done
254
255 // resu is set to zero in the other domains :
256
257 if (l_min > 0) {
258 resu.annule(0, l_min-1) ;
259 }
260
261 if (l_min + nzet < nz) {
262 resu.annule(l_min + nzet, nz - 1) ;
263 }
264}
265
266
267
268void Eos::calcule(const Scalar& ent, int nzet, int l_min,
269 double (Eos::*fait)(double, const Param*) const, const Param* par, Scalar& resu) const {
270
271 assert(ent.get_etat() != ETATNONDEF) ;
272
273 const Map* mp = &(ent.get_mp()) ; // Mapping
274
275
276 if (ent.get_etat() == ETATZERO) {
277 resu.set_etat_zero() ;
278 return ;
279 }
280
281 assert(ent.get_etat() == ETATQCQ) ;
282 const Valeur& vent = ent.get_spectral_va() ;
283 vent.coef_i() ; // the values in the configuration space are required
284
285 const Mg3d* mg = mp->get_mg() ; // Multi-grid
286
287 int nz = mg->get_nzone() ; // total number of domains
288
289 // Preparations for a point by point computation:
290 resu.set_etat_qcq() ;
291 Valeur& vresu = resu.set_spectral_va() ;
292 vresu.set_etat_c_qcq() ;
293 vresu.c->set_etat_qcq() ;
294
295 // Loop on domains where the computation has to be done :
296 for (int l = l_min; l< l_min + nzet; l++) {
297
298 assert(l>=0) ;
299 assert(l<nz) ;
300
301 Tbl* tent = vent.c->t[l] ;
302 Tbl* tresu = vresu.c->t[l] ;
303
304 if (tent->get_etat() == ETATZERO) {
305 tresu->set_etat_zero() ;
306 }
307 else {
308 assert( tent->get_etat() == ETATQCQ ) ;
309 tresu->set_etat_qcq() ;
310
311 for (int i=0; i<tent->get_taille(); i++) {
312
313 tresu->t[i] = (this->*fait)( tent->t[i], par ) ;
314 }
315
316 } // End of the case where ent != 0 in the considered domain
317
318 } // End of the loop on domains where the computation had to be done
319
320 // resu is set to zero in the other domains :
321
322 if (l_min > 0) {
323 resu.annule(0, l_min-1) ;
324 }
325
326 if (l_min + nzet < nz) {
327 resu.annule(l_min + nzet, nz - 1) ;
328 }
329}
330 //---------------------------------//
331 // Public functions //
332 //---------------------------------//
333
334
335// Baryon density from enthalpy
336//------------------------------
337
338Cmp Eos::nbar_ent(const Cmp& ent, int nzet, int l_min, const Param* par) const {
339
340 Cmp resu(ent.get_mp()) ;
341
342 calcule(ent, nzet, l_min, &Eos::nbar_ent_p, par, resu) ;
343
344 return resu ;
345
346}
347
348Scalar Eos::nbar_ent(const Scalar& ent, int nzet, int l_min, const Param* par) const {
349
350 Scalar resu(ent.get_mp()) ;
351
352 calcule(ent, nzet, l_min, &Eos::nbar_ent_p, par, resu) ;
353
354 return resu ;
355
356}
357
358
359
360// Energy density from enthalpy
361//------------------------------
362
363Cmp Eos::ener_ent(const Cmp& ent, int nzet, int l_min, const Param* par) const {
364
365 Cmp resu(ent.get_mp()) ;
366
367 calcule(ent, nzet, l_min, &Eos::ener_ent_p, par, resu) ;
368
369 return resu ;
370
371}
372
373Scalar Eos::ener_ent(const Scalar& ent, int nzet, int l_min, const Param* par) const {
374
375 Scalar resu(ent.get_mp()) ;
376
377 calcule(ent, nzet, l_min, &Eos::ener_ent_p, par, resu) ;
378
379 return resu ;
380
381}
382// Pressure from enthalpy
383//-----------------------
384
385Cmp Eos::press_ent(const Cmp& ent, int nzet, int l_min, const Param* par) const {
386
387 Cmp resu(ent.get_mp()) ;
388
389 calcule(ent, nzet, l_min, &Eos::press_ent_p, par, resu) ;
390
391 return resu ;
392
393}
394
395Scalar Eos::press_ent(const Scalar& ent, int nzet, int l_min, const Param* par) const {
396
397 Scalar resu(ent.get_mp()) ;
398
399 calcule(ent, nzet, l_min, &Eos::press_ent_p, par, resu) ;
400
401 return resu ;
402
403}
404// Derivative of baryon density from enthalpy
405//-------------------------------------------
406
407Cmp Eos::der_nbar_ent(const Cmp& ent, int nzet, int l_min, const Param* par) const {
408
409 Cmp resu(ent.get_mp()) ;
410
411 calcule(ent, nzet, l_min, &Eos::der_nbar_ent_p, par, resu) ;
412
413 return resu ;
414
415}
416
417Scalar Eos::der_nbar_ent(const Scalar& ent, int nzet, int l_min, const Param* par) const {
418
419 Scalar resu(ent.get_mp()) ;
420
421 calcule(ent, nzet, l_min, &Eos::der_nbar_ent_p, par, resu) ;
422
423 return resu ;
424
425}
426
427// Derivative of energy density from enthalpy
428//-------------------------------------------
429
430Cmp Eos::der_ener_ent(const Cmp& ent, int nzet, int l_min, const Param* par) const {
431
432 Cmp resu(ent.get_mp()) ;
433
434 calcule(ent, nzet, l_min, &Eos::der_ener_ent_p, par, resu) ;
435
436 return resu ;
437
438}
439
440Scalar Eos::der_ener_ent(const Scalar& ent, int nzet, int l_min, const Param* par) const {
441
442 Scalar resu(ent.get_mp()) ;
443
444 calcule(ent, nzet, l_min, &Eos::der_ener_ent_p, par, resu) ;
445
446 return resu ;
447
448}
449// Derivative of pressure from enthalpy
450//-------------------------------------------
451
452Cmp Eos::der_press_ent(const Cmp& ent, int nzet, int l_min, const Param* par) const {
453
454 Cmp resu(ent.get_mp()) ;
455
456 calcule(ent, nzet, l_min, &Eos::der_press_ent_p, par, resu) ;
457
458 return resu ;
459
460}
461
462Scalar Eos::der_press_ent(const Scalar& ent, int nzet, int l_min, const Param* par) const {
463
464 Scalar resu(ent.get_mp()) ;
465
466 calcule(ent, nzet, l_min, &Eos::der_press_ent_p, par, resu) ;
467
468 return resu ;
469
470}
471}
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition cmp.h:446
int get_etat() const
Returns the logical state.
Definition cmp.h:899
Valeur va
The numerical value of the Cmp
Definition cmp.h:464
const Map * get_mp() const
Returns the mapping.
Definition cmp.h:901
Equation of state base class.
Definition eos.h:190
virtual int identify() const =0
Returns a number to identify the sub-classe of Eos the object belongs to.
void calcule(const Cmp &thermo, int nzet, int l_min, double(Eos::*fait)(double, const Param *) const, const Param *par, Cmp &resu) const
General computational method for Cmp 's.
Definition eos.C:203
virtual ~Eos()
Destructor.
Definition eos.C:152
virtual void sauve(FILE *) const
Save in a file.
Definition eos.C:179
const char * get_name() const
Returns the EOS name.
Definition eos.C:169
Eos()
Standard constructor.
Definition eos.C:108
char name[100]
EOS name.
Definition eos.h:196
void set_name(const char *name_i)
Sets the EOS name.
Definition eos.C:163
Base class for coordinate mappings.
Definition map.h:670
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition map.h:765
Multi-domain grid.
Definition grilles.h:273
int get_nzone() const
Returns the number of domains.
Definition grilles.h:448
Tbl ** t
Array (size nzone ) of pointers on the Tbl 's.
Definition mtbl.h:132
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition mtbl.C:299
Parameter storage.
Definition param.h:125
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
const Valeur & get_spectral_va() const
Returns va (read only version)
Definition scalar.h:601
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
Definition scalar.h:554
Basic array class.
Definition tbl.h:161
int get_etat() const
Gives the logical state.
Definition tbl.h:394
void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition tbl.C:347
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition tbl.C:361
int get_taille() const
Gives the total size (ie dim.taille)
Definition tbl.h:397
double * t
The array of double.
Definition tbl.h:173
Values and coefficients of a (real-value) function.
Definition valeur.h:287
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
Mtbl * c
Values of the function at the points of the multi-grid
Definition valeur.h:299
void coef_i() const
Computes the physical value of *this.
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
const Map & get_mp() const
Returns the mapping.
Definition tensor.h:861
Lorene prototypes.
Definition app_hor.h:64