29 char eos_bifluid_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Eos/eos_bifluid.C,v 1.22 2015/06/12 12:38:24 j_novak Exp $" ;
136 #include "eos_bifluid.h" 138 #include "utilitaires.h" 149 name(
"EoS bi-fluid"), m_1(1), m_2(1)
168 char dummy [MAX_EOSNAME];
169 if (fread(dummy,
sizeof(
char),MAX_EOSNAME, fich) == 0)
170 cerr <<
"Reading problem in " << __FILE__ << endl ;
183 char* char_name =
const_cast<char*
>(
"name");
184 char* char_m1 =
const_cast<char*
>(
"m_1") ;
185 char* char_m2 =
const_cast<char*
>(
"m_2") ;
196 cerr <<
"ERROR: Eos_bifluid(const char*): could not read either of \ 197 the variables 'name', 'm_1, or 'm_2' from file " << fname << endl;
213 fich.getline(blabla, 80) ;
215 fich >>
m_1 ; fich.getline(blabla, 80) ;
216 fich >>
m_2 ; fich.getline(blabla, 80) ;
247 assert(
name.size() < MAX_EOSNAME) ;
248 char dummy [MAX_EOSNAME];
251 fwrite_be(&ident,
sizeof(
int), 1, fich) ;
253 strncpy (dummy,
name.c_str(), MAX_EOSNAME);
254 dummy[MAX_EOSNAME-1] = 0;
255 if (fwrite(dummy,
sizeof(
char), MAX_EOSNAME, fich ) == 0)
256 cerr <<
"Writing problem in " << __FILE__ << endl ;
266 ost <<
" Mean particle 1 mass : " << eqetat.
get_m1() <<
" m_B" << endl ;
267 ost <<
" Mean particle 2 mass : " << eqetat.
get_m2() <<
" m_B" << endl ;
282 const Cmp& delta2,
Cmp& nbar1,
Cmp& nbar2,
284 int nzet,
int l_min)
const {
286 assert(ent1.
get_etat() != ETATNONDEF) ;
287 assert(ent2.
get_etat() != ETATNONDEF) ;
288 assert(delta2.
get_etat() != ETATNONDEF) ;
291 assert(mp == ent2.
get_mp()) ;
292 assert(mp == delta2.
get_mp()) ;
293 assert(mp == ener.
get_mp()) ;
307 const Mg3d* mg = mp->get_mg() ;
314 nbar1.
annule(0, l_min-1) ;
315 nbar2.
annule(0, l_min-1) ;
317 press.
annule(0, l_min-1) ;
320 if (l_min + nzet < nz) {
321 nbar1.
annule(l_min + nzet, nz - 1) ;
322 nbar2.
annule(l_min + nzet, nz - 1) ;
323 ener.
annule(l_min + nzet, nz - 1) ;
324 press.
annule(l_min + nzet, nz - 1) ;
327 int np0 = mp->get_mg()->get_np(0) ;
328 int nt0 = mp->get_mg()->get_nt(0) ;
329 for (
int l=l_min; l<l_min+nzet; l++) {
330 assert(mp->get_mg()->get_np(l) == np0) ;
331 assert(mp->get_mg()->get_nt(l) == nt0) ;
337 bool slow_rot_style =
false;
342 if (this_eos -> get_typeos() == 5)
344 slow_rot_style =
true;
345 cout <<
"DEBUG: using EOS-inversion slow-rot-style!! " << endl;
351 for (
int k=0; k<np0; k++) {
352 for (
int j=0; j<nt0; j++) {
354 bool inside1 = true ;
355 bool inside2 = true ;
356 for (
int l=l_min; l< l_min + nzet; l++) {
357 for (
int i=0; i<mp->get_mg()->get_nr(l); i++) {
358 double xx, xx1, xx2, n1, n2 ;
359 xx1 = ent1(l,k,j,i) ;
360 xx2 = ent2(l,k,j,i) ;
361 xx = delta2(l,k,j,i) ;
363 inside = (!
nbar_ent_p(xx1, xx2, xx, n1, n2)) ;
366 inside1 = (n1 > 0.) ;
368 inside2 = (n2 > 0.) ;
374 n1 = (n1 > 0) ? n1: 0;
375 n2 = (n2 > 0) ? n2: 0;
379 nbar1.
set(l,k,j,i) = n1 ;
380 nbar2.
set(l,k,j,i) = n2 ;
385 inside1 = (n1 > 0.) ;
389 inside2 = (n2 > 0.) ;
391 if (!inside1) n1 = 0. ;
392 if (!inside2) n2 = 0. ;
393 nbar1.
set(l,k,j,i) = n1 ;
394 nbar2.
set(l,k,j,i) = n2 ;
412 Cmp& nbar1,
Cmp& nbar2,
int nzet,
int l_min)
const {
414 assert(ent1.
get_etat() != ETATNONDEF) ;
415 assert(ent2.
get_etat() != ETATNONDEF) ;
416 assert(delta2.
get_etat() != ETATNONDEF) ;
419 assert(mp == ent2.
get_mp()) ;
420 assert(mp == delta2.
get_mp()) ;
421 assert(mp == nbar1.
get_mp()) ;
431 const Mg3d* mg = mp->get_mg() ;
438 nbar1.
annule(0, l_min-1) ;
439 nbar2.
annule(0, l_min-1) ;
442 if (l_min + nzet < nz) {
443 nbar1.
annule(l_min + nzet, nz - 1) ;
444 nbar2.
annule(l_min + nzet, nz - 1) ;
447 int np0 = mp->get_mg()->get_np(0) ;
448 int nt0 = mp->get_mg()->get_nt(0) ;
449 for (
int l=l_min; l<l_min+nzet; l++) {
450 assert(mp->get_mg()->get_np(l) == np0) ;
451 assert(mp->get_mg()->get_nt(l) == nt0) ;
454 for (
int k=0; k<np0; k++) {
455 for (
int j=0; j<nt0; j++) {
457 bool inside1 = true ;
458 bool inside2 = true ;
459 for (
int l=l_min; l< l_min + nzet; l++) {
460 for (
int i=0; i<mp->get_mg()->get_nr(l); i++) {
461 double xx, xx1, xx2, n1, n2 ;
462 xx1 = ent1(l,k,j,i) ;
463 xx2 = ent2(l,k,j,i) ;
464 xx = delta2(l,k,j,i) ;
466 inside = (!
nbar_ent_p(xx1, xx2, xx, n1, n2)) ;
467 inside1 = ((n1 > 0.)&&(xx1>0.)) ;
468 inside2 = ((n2 > 0.)&&(xx2>0.)) ;
471 nbar1.
set(l,k,j,i) = n1 ;
472 nbar2.
set(l,k,j,i) = n2 ;
477 inside1 = (n1 > 0.) ;
479 if (!inside1) n1 = 0. ;
482 inside2 = (n2 > 0.) ;
484 if (!inside2) n2 = 0. ;
485 nbar1.
set(l,k,j,i) = n1 ;
486 nbar2.
set(l,k,j,i) = n2 ;
500 int nzet,
int l_min)
const {
503 assert(ent1.
get_etat() != ETATNONDEF) ;
504 assert(ent2.
get_etat() != ETATNONDEF) ;
505 assert(delta2.
get_etat() != ETATNONDEF) ;
508 assert(mp == ent2.
get_mp()) ;
509 assert(mp == delta2.
get_mp()) ;
512 ener.set_etat_zero() ;
516 ener.allocate_all() ;
518 const Mg3d* mg = mp->get_mg() ;
525 ener.annule(0, l_min-1) ;
527 if (l_min + nzet < nz)
528 ener.annule(l_min + nzet, nz - 1) ;
530 int np0 = mp->get_mg()->get_np(0) ;
531 int nt0 = mp->get_mg()->get_nt(0) ;
532 for (
int l=l_min; l<l_min+nzet; l++) {
533 assert(mp->get_mg()->get_np(l) == np0) ;
534 assert(mp->get_mg()->get_nt(l) == nt0) ;
537 for (
int k=0; k<np0; k++) {
538 for (
int j=0; j<nt0; j++) {
540 bool inside1 = true ;
541 bool inside2 = true ;
542 for (
int l=l_min; l< l_min + nzet; l++) {
543 for (
int i=0; i<mp->get_mg()->get_nr(l); i++) {
544 double xx, xx1, xx2, n1, n2 ;
545 xx1 = ent1(l,k,j,i) ;
546 xx2 = ent2(l,k,j,i) ;
547 xx = delta2(l,k,j,i) ;
549 inside = (!
nbar_ent_p(xx1, xx2, xx, n1, n2)) ;
550 inside1 = ((n1 > 0.)&&(xx1>0.)) ;
551 inside2 = ((n2 > 0.)&&(xx2>0.)) ;
559 inside1 = (n1 > 0.) ;
561 if (!inside1) n1 = 0. ;
564 inside2 = (n2 > 0.) ;
566 if (!inside2) n2 = 0. ;
580 int nzet,
int l_min )
const {
583 assert(ent1.
get_etat() != ETATNONDEF) ;
584 assert(ent2.
get_etat() != ETATNONDEF) ;
585 assert(delta2.
get_etat() != ETATNONDEF) ;
588 assert(mp == ent2.
get_mp()) ;
591 press.set_etat_zero() ;
594 press.allocate_all() ;
596 const Mg3d* mg = mp->get_mg() ;
603 press.annule(0, l_min-1) ;
605 if (l_min + nzet < nz)
606 press.annule(l_min + nzet, nz - 1) ;
608 int np0 = mp->get_mg()->get_np(0) ;
609 int nt0 = mp->get_mg()->get_nt(0) ;
610 for (
int l=l_min; l<l_min+nzet; l++) {
611 assert(mp->get_mg()->get_np(l) == np0) ;
612 assert(mp->get_mg()->get_nt(l) == nt0) ;
615 for (
int k=0; k<np0; k++) {
616 for (
int j=0; j<nt0; j++) {
618 bool inside1 = true ;
619 bool inside2 = true ;
620 for (
int l=l_min; l< l_min + nzet; l++) {
621 for (
int i=0; i<mp->get_mg()->get_nr(l); i++) {
622 double xx, xx1, xx2, n1, n2 ;
623 xx1 = ent1(l,k,j,i) ;
624 xx2 = ent2(l,k,j,i) ;
625 xx = delta2(l,k,j,i) ;
627 inside = (!
nbar_ent_p(xx1, xx2, xx, n1, n2)) ;
628 inside1 = ((n1 > 0.)&&(xx1>0.)) ;
629 inside2 = ((n2 > 0.)&&(xx2>0.)) ;
636 inside1 = (n1 > 0.) ;
638 if (!inside1) n1 = 0. ;
641 inside2 = (n2 > 0.) ;
643 if (!inside2) n2 = 0. ;
657 x2,
int nzet,
int l_min,
double 658 (
Eos_bifluid::*fait)(
double,
double,
double)
const,
661 assert(nbar1.
get_etat() != ETATNONDEF) ;
662 assert(nbar2.
get_etat() != ETATNONDEF) ;
663 assert(x2.
get_etat() != ETATNONDEF) ;
666 assert(mp == nbar2.
get_mp()) ;
674 bool nb1 = nbar1.
get_etat() == ETATQCQ ;
675 bool nb2 = nbar2.
get_etat() == ETATQCQ ;
676 bool bx2 = x2.
get_etat() == ETATQCQ ;
677 const Valeur* vnbar1 = 0x0 ;
678 const Valeur* vnbar2 = 0x0 ;
680 if (nb1) { vnbar1 = &nbar1.
va ;
682 if (nb2) { vnbar2 = &nbar2.
va ;
684 if (bx2) {vx2 = & x2.
va ;
687 const Mg3d* mg = mp->get_mg() ;
698 for (
int l = l_min; l< l_min + nzet; l++) {
707 if (nb1) tnbar1 = vnbar1->
c->
t[l] ;
708 if (nb2) tnbar2 = vnbar2->
c->
t[l] ;
709 if (bx2) tx2 = vx2->
c->
t[l] ;
710 Tbl* tresu = vresu.
c->
t[l] ;
713 if (nb1) nb1b = tnbar1->
get_etat() == ETATQCQ ;
715 if (nb2) nb2b = tnbar2->
get_etat() == ETATQCQ ;
717 if (bx2) bx2b = tx2->
get_etat() == ETATQCQ ;
724 n1 = nb1b ? tnbar1->
t[i] : 0 ;
725 n2 = nb2b ? tnbar2->
t[i] : 0 ;
726 xx = bx2b ? tx2->
t[i] : 0 ;
727 tresu->
t[i] = (this->*fait)(n1, n2, xx ) ;
740 if (l_min + nzet < nz) {
741 resu.
annule(l_min + nzet, nz - 1) ;
746 delta2,
int nzet,
int l_min)
const {
757 delta2,
int nzet,
int l_min)
const {
768 delta2,
int nzet,
int l_min)
const {
784 x2,
int nzet,
int l_min,
double 785 (
Eos_bifluid::*fait)(
double,
double,
double)
const,
788 assert(ent1.
get_etat() != ETATNONDEF) ;
789 assert(ent2.
get_etat() != ETATNONDEF) ;
790 assert(x2.
get_etat() != ETATNONDEF) ;
793 assert(mp == ent2.
get_mp()) ;
801 bool le1 = ent1.
get_etat() == ETATQCQ ;
802 bool le2 = ent2.
get_etat() == ETATQCQ ;
803 bool bx2 = x2.
get_etat() == ETATQCQ ;
804 const Valeur* vent1 = 0x0 ;
805 const Valeur* vent2 = 0x0 ;
807 if (le1) { vent1 = &ent1.
va ;
809 if (le2) { vent2 = &ent2.
va ;
811 if (bx2) {vx2 = & x2.
va ;
814 const Mg3d* mg = mp->get_mg() ;
825 for (
int l = l_min; l< l_min + nzet; l++) {
834 if (le1) tent1 = vent1->
c->
t[l] ;
835 if (le2) tent2 = vent2->
c->
t[l] ;
836 if (bx2) tx2 = vx2->
c->
t[l] ;
837 Tbl* tresu = vresu.
c->
t[l] ;
840 if (le1) le1b = tent1->
get_etat() == ETATQCQ ;
842 if (le2) le2b = tent2->
get_etat() == ETATQCQ ;
844 if (bx2) bx2b = tx2->
get_etat() == ETATQCQ ;
851 H1 = le1b ? tent1->
t[i] : 0 ;
852 H2 = le2b ? tent2->
t[i] : 0 ;
853 xx = bx2b ? tx2->
t[i] : 0 ;
856 tresu->
t[i] = (this->*fait)(xx, H1, H2) ;
867 if (l_min + nzet < nz) {
868 resu.
annule(l_min + nzet, nz - 1) ;
874 delta2,
int nzet,
int l_min)
const {
885 delta2,
int nzet,
int l_min)
const {
896 delta2,
int nzet,
int l_min)
const {
Cmp get_Kpp(const Cmp &nbar1, const Cmp &nbar2, const Cmp &x2, int nzet, int l_min=0) const
Computes the derivatives of the energy/(baryonic density 2) .
virtual double get_K12(const double n1, const double n2, const double x) const =0
Computes the derivative of the energy with respect to .
const Map * get_mp() const
Returns the mapping.
double get_m2() const
Return the individual particule mass .
virtual void sauve(FILE *) const
Save in a file.
virtual double get_K22(const double n1, const double n2, const double x) const =0
Computes the derivative of the energy/(baryonic density 2) .
double m_1
Individual particle mass [unit: ].
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Cmp press_ent(const Cmp &ent1, const Cmp &ent2, const Cmp &delta2, int nzet, int l_min=0) const
Computes the pressure from the log-enthalpy fields and the relative velocity.
virtual double nbar_ent_p1(const double ent1) const =0
Computes baryon density out of the log-enthalpy asuming that only fluid 1 is present (virtual functio...
void annule(int l)
Sets the Cmp to zero in a given domain.
Cmp get_Knn_ent(const Cmp &ent1, const Cmp &ent2, const Cmp &x2, int nzet, int l_min=0) const
Computes the derivatives of the energy/(baryonic density 1) .
virtual int identify() const =0
Returns a number to identify the sub-classe of Eos_bifluid the object belongs to. ...
Cmp ener_ent(const Cmp &ent1, const Cmp &ent2, const Cmp &delta2, int nzet, int l_min=0) const
Computes the total energy density from the log-enthalpy fields and the relative velocity.
int get_etat() const
Returns the logical state.
void coef_i() const
Computes the physical value of *this.
Base class for coordinate mappings.
Values and coefficients of a (real-value) function.
int get_etat() const
Gives the logical state.
void operator=(const Eos_bifluid &)
Assignment to another Eos_bifluid.
2-fluids equation of state base class.
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Cmp get_Knn(const Cmp &nbar1, const Cmp &nbar2, const Cmp &x2, int nzet, int l_min=0) const
Computes the derivatives of the energy/(baryonic density 1) .
virtual double get_K11(const double n1, const double n2, const double x) const =0
Computes the derivative of the energy with respect to (baryonic density 1) .
double m_2
Individual particle mass [unit: ].
Cmp get_Kpp_ent(const Cmp &ent1, const Cmp &ent2, const Cmp &x2, int nzet, int l_min=0) const
Computes the derivatives of the energy/(baryonic density 2) .
void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Eos_bifluid()
Standard constructor.
Cmp get_Knp_ent(const Cmp &ent1, const Cmp &ent2, const Cmp &x2, int nzet, int l_min=0) const
Computes the derivatives of the energy with respect to .
virtual double ener_nbar_p(const double nbar1, const double nbar2, const double delta2) const =0
Computes the total energy density from the baryonic densities and the relative velocity.
void nbar_ent(const Cmp &ent1, const Cmp &ent2, const Cmp &delta2, Cmp &nbar1, Cmp &nbar2, int nzet, int l_min=0) const
Computes both baryon density fields from the log-enthalpy fields and the relative velocity...
double * t
The array of double.
virtual bool nbar_ent_p(const double ent1, const double ent2, const double delta2, double &nbar1, double &nbar2) const =0
Computes both baryon densities from the log-enthalpies (virtual function implemented in the derived c...
Mtbl * c
Values of the function at the points of the multi-grid.
int get_nzone() const
Returns the number of domains.
Analytic equation of state for two fluids (Newtonian case).
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.
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
friend ostream & operator<<(ostream &, const Eos_bifluid &)
Display.
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
virtual double nbar_ent_p2(const double ent2) const =0
Computes baryon density out of the log-enthalpy assuming that only fluid 2 is present (virtual functi...
void calcule(const Cmp &nbar1, const Cmp &nbar2, const Cmp &x2, int nzet, int l_min, double(Eos_bifluid::*fait)(double, double, double) const, Cmp &resu) const
General computational method for Cmp 's ( 's).
virtual double press_nbar_p(const double nbar1, const double nbar2, const double delta2) const =0
Computes the pressure from the baryonic densities and the relative velocity.
int fread_be(int *aa, int size, int nb, FILE *fich)
Reads integer(s) from a binary file according to the big endian convention.
void allocate_all()
Sets the logical state to ETATQCQ (ordinary state) and performs the memory allocation of all the elem...
Tbl & set(int l)
Read/write of the value in a given domain.
virtual void calcule_tout(const Cmp &ent1, const Cmp &ent2, const Cmp &delta2, Cmp &nbar1, Cmp &nbar2, Cmp &ener, Cmp &press, int nzet, int l_min=0) const
General computational method for Cmp 's, it computes both baryon densities, energy and pressure profi...
virtual ~Eos_bifluid()
Destructor.
string get_name() const
Returns the EOS name.
int read_variable(const char *fname, const char *var_name, char *fmt, void *varp)
Reads a variable from file.
int get_taille() const
Gives the total size (ie dim.taille)
void set_etat_c_qcq()
Sets the logical state to ETATQCQ (ordinary state) for values in the configuration space (Mtbl c )...
Cmp get_Knp(const Cmp &nbar1, const Cmp &nbar2, const Cmp &x2, int nzet, int l_min=0) const
Computes the derivatives of the energy with respect to .
Tbl ** t
Array (size nzone ) of pointers on the Tbl 's.
Valeur va
The numerical value of the Cmp.
double get_m1() const
Return the individual particule mass .
void calcule_ent(const Cmp &ent1, const Cmp &ent2, const Cmp &x2, int nzet, int l_min, double(Eos_bifluid::*fait)(double, double, double) const, Cmp &resu) const
General computational method for Cmp 's ( 's).