valeur_math.C
/*
* Mathematical functions for class Valeur
*
*/
/*
* Copyright (c) 1999-2000 Jean-Alain Marck
* Copyright (c) 1999-2001 Eric Gourgoulhon
* Copyright (c) 1999-2001 Philippe Grandclement
*
* This file is part of LORENE.
*
* LORENE is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation; either version 2 of the License, or
* (at your option) any later version.
*
* LORENE is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with LORENE; if not, write to the Free Software
* Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
*
*/
char valeur_math_C[] = "$Header: /cvsroot/Lorene/C++/Source/Valeur/valeur_math.C,v 1.5 2014/10/13 08:53:50 j_novak Exp $" ;
/*
* $Id: valeur_math.C,v 1.5 2014/10/13 08:53:50 j_novak Exp $
* $Log: valeur_math.C,v $
* Revision 1.5 2014/10/13 08:53:50 j_novak
* Lorene classes and functions now belong to the namespace Lorene.
*
* Revision 1.4 2014/10/06 15:13:23 j_novak
* Modified #include directives to use c++ syntax.
*
* Revision 1.3 2012/01/17 10:39:27 j_penner
* added a Heaviside function
*
* Revision 1.2 2002/10/16 14:37:15 j_novak
* Reorganization of #include instructions of standard C++, in order to
* use experimental version 3 of gcc.
*
* Revision 1.1.1.1 2001/11/20 15:19:27 e_gourgoulhon
* LORENE
*
* Revision 2.4 1999/12/02 17:57:42 phil
* *** empty log message ***
*
* Revision 2.3 1999/11/09 16:18:02 phil
* ajout de racine_cubique
*
* Revision 2.2 1999/10/29 15:14:50 eric
* Ajout de fonctions mathematiques (abs, norme, etc...).
*
* Revision 2.1 1999/03/01 15:11:03 eric
* *** empty log message ***
*
* Revision 2.0 1999/02/24 15:40:32 hyc
* *** empty log message ***
*
*
* $Header: /cvsroot/Lorene/C++/Source/Valeur/valeur_math.C,v 1.5 2014/10/13 08:53:50 j_novak Exp $
*
*/
// Fichiers include
// ----------------
#include <cmath>
#include <cstdlib>
#include "valeur.h"
//-------//
// Sinus //
//-------//
namespace Lorene {
Valeur sin(const Valeur& ti)
{
// Protection
assert(ti.get_etat() != ETATNONDEF) ;
// Cas ETATZERO
if (ti.get_etat() == ETATZERO) {
return ti ;
}
// Cas general
assert(ti.get_etat() == ETATQCQ) ; // sinon...
if (ti.c == 0x0) { // Il faut la valeur physique
ti.coef_i() ;
}
Valeur to(ti.get_mg()) ; // Valeur resultat
to.set_etat_c_qcq() ;
*(to.c) = sin( *(ti.c) ) ;
return to ;
}
//---------//
// Cosinus //
//---------//
Valeur cos(const Valeur& ti)
{
// Protection
assert(ti.get_etat() != ETATNONDEF) ;
Valeur to(ti.get_mg()) ; // Valeur resultat
to.set_etat_c_qcq() ;
// Cas ETATZERO
if (ti.get_etat() == ETATZERO) {
*(to.c) = 1. ;
return to ;
}
// Cas general
assert(ti.get_etat() == ETATQCQ) ; // sinon...
if (ti.c == 0x0) { // Il faut la valeur physique
ti.coef_i() ;
}
*(to.c) = cos( *(ti.c) ) ;
return to ;
}
//----------//
// Tangente //
//----------//
Valeur tan(const Valeur& ti)
{
// Protection
assert(ti.get_etat() != ETATNONDEF) ;
// Cas ETATZERO
if (ti.get_etat() == ETATZERO) {
return ti ;
}
// Cas general
assert(ti.get_etat() == ETATQCQ) ; // sinon...
if (ti.c == 0x0) { // Il faut la valeur physique
ti.coef_i() ;
}
Valeur to(ti.get_mg()) ; // Valeur resultat
to.set_etat_c_qcq() ;
*(to.c) = tan( *(ti.c) ) ;
return to ;
}
//----------//
// ArcSinus //
//----------//
Valeur asin(const Valeur& ti)
{
// Protection
assert(ti.get_etat() != ETATNONDEF) ;
// Cas ETATZERO
if (ti.get_etat() == ETATZERO) {
return ti ;
}
// Cas general
assert(ti.get_etat() == ETATQCQ) ; // sinon...
if (ti.c == 0x0) { // Il faut la valeur physique
ti.coef_i() ;
}
Valeur to(ti.get_mg()) ; // Valeur resultat
to.set_etat_c_qcq() ;
*(to.c) = asin( *(ti.c) ) ;
return to ;
}
//------------//
// ArcCosinus //
//------------//
Valeur acos(const Valeur& ti)
{
// Protection
assert(ti.get_etat() != ETATNONDEF) ;
Valeur to(ti.get_mg()) ; // Valeur resultat
to.set_etat_c_qcq() ;
// Cas ETATZERO
if (ti.get_etat() == ETATZERO) {
*(to.c) = M_PI * .5 ;
return to ;
}
// Cas general
assert(ti.get_etat() == ETATQCQ) ; // sinon...
if (ti.c == 0x0) { // Il faut la valeur physique
ti.coef_i() ;
}
*(to.c) = acos( *(ti.c) ) ;
return to ;
}
//-------------//
// ArcTangente //
//-------------//
Valeur atan(const Valeur& ti)
{
// Protection
assert(ti.get_etat() != ETATNONDEF) ;
// Cas ETATZERO
if (ti.get_etat() == ETATZERO) {
return ti ;
}
// Cas general
assert(ti.get_etat() == ETATQCQ) ; // sinon...
if (ti.c == 0x0) { // Il faut la valeur physique
ti.coef_i() ;
}
Valeur to(ti.get_mg()) ; // Valeur resultat
to.set_etat_c_qcq() ;
*(to.c) = atan( *(ti.c) ) ;
return to ;
}
//------//
// Sqrt //
//------//
Valeur sqrt(const Valeur& ti)
{
// Protection
assert(ti.get_etat() != ETATNONDEF) ;
// Cas ETATZERO
if (ti.get_etat() == ETATZERO) {
return ti ;
}
// Cas general
assert(ti.get_etat() == ETATQCQ) ; // sinon...
if (ti.c == 0x0) { // Il faut la valeur physique
ti.coef_i() ;
}
Valeur to(ti.get_mg()) ; // Valeur resultat
to.set_etat_c_qcq() ;
*(to.c) = sqrt( *(ti.c) ) ;
return to ;
}
//---------------//
// Exponantielle //
//---------------//
Valeur exp(const Valeur& ti)
{
// Protection
assert(ti.get_etat() != ETATNONDEF) ;
Valeur to(ti.get_mg()) ; // Valeur resultat
to.set_etat_c_qcq() ;
// Cas ETATZERO
if (ti.get_etat() == ETATZERO) {
*(to.c) = 1. ;
return to ;
}
// Cas general
assert(ti.get_etat() == ETATQCQ) ; // sinon...
if (ti.c == 0x0) { // Il faut la valeur physique
ti.coef_i() ;
}
*(to.c) = exp( *(ti.c) ) ;
return to ;
}
//--------------------//
// Heaviside Function //
//--------------------//
Valeur Heaviside(const Valeur& ti)
{
// Protection
assert(ti.get_etat() != ETATNONDEF) ;
Valeur to(ti.get_mg()) ; // Valeur resultat
to.set_etat_c_qcq() ;
// Cas ETATZERO
if (ti.get_etat() == ETATZERO) {
*(to.c) = 0. ;
return to ;
}
// Cas general
assert(ti.get_etat() == ETATQCQ) ; // otherwise
if (ti.c == 0x0) { // Use the physical value << What?? (used Google translate)
ti.coef_i() ;
}
*(to.c) = Heaviside( *(ti.c) ) ;
return to ;
}
//-------------//
// Log naturel //
//-------------//
Valeur log(const Valeur& ti)
{
// Protection
assert(ti.get_etat() != ETATNONDEF) ;
// Cas ETATZERO
if (ti.get_etat() == ETATZERO) {
cout << "Valeur log: log(ETATZERO) !" << endl ;
abort () ;
}
// Cas general
assert(ti.get_etat() == ETATQCQ) ; // sinon...
if (ti.c == 0x0) { // Il faut la valeur physique
ti.coef_i() ;
}
Valeur to(ti.get_mg()) ; // Valeur resultat
to.set_etat_c_qcq() ;
*(to.c) = log( *(ti.c) ) ;
return to ;
}
//-------------//
// Log decimal //
//-------------//
Valeur log10(const Valeur& ti)
{
// Protection
assert(ti.get_etat() != ETATNONDEF) ;
// Cas ETATZERO
if (ti.get_etat() == ETATZERO) {
cout << "Valeur log10: log10(ETATZERO) !" << endl ;
abort () ;
}
// Cas general
assert(ti.get_etat() == ETATQCQ) ; // sinon...
if (ti.c == 0x0) { // Il faut la valeur physique
ti.coef_i() ;
}
Valeur to(ti.get_mg()) ; // Valeur resultat
to.set_etat_c_qcq() ;
*(to.c) = log10( *(ti.c) ) ;
return to ;
}
//--------------//
// Power entier //
//--------------//
Valeur pow(const Valeur& ti, int n)
{
// Protection
assert(ti.get_etat() != ETATNONDEF) ;
// Cas ETATZERO
if (ti.get_etat() == ETATZERO) {
if (n > 0) {
return ti ;
}
else {
cout << "Valeur pow: ETATZERO^n with n<=0 ! "<< endl ;
abort () ;
}
}
// Cas general
assert(ti.get_etat() == ETATQCQ) ; // sinon...
if (ti.c == 0x0) { // Il faut la valeur physique
ti.coef_i() ;
}
Valeur to(ti.get_mg()) ; // Valeur resultat
to.set_etat_c_qcq() ;
double x = n ;
*(to.c) = pow( *(ti.c), x ) ;
return to ;
}
//--------------//
// Power double //
//--------------//
Valeur pow(const Valeur& ti, double x)
{
// Protection
assert(ti.get_etat() != ETATNONDEF) ;
// Cas ETATZERO
if (ti.get_etat() == ETATZERO) {
if (x > 0) {
return ti ;
}
else {
cout << "Valeur pow: ETATZERO^x with x<=0 !" << endl ;
abort () ;
}
}
// Cas general
assert(ti.get_etat() == ETATQCQ) ; // sinon...
if (ti.c == 0x0) { // Il faut la valeur physique
ti.coef_i() ;
}
Valeur to(ti.get_mg()) ; // Valeur resultat
to.set_etat_c_qcq() ;
*(to.c) = pow( *(ti.c), x ) ;
return to ;
}
//----------------//
// Absolute value //
//----------------//
Valeur abs(const Valeur& vi)
{
// Protection
assert(vi.get_etat() != ETATNONDEF) ;
// Cas ETATZERO
if (vi.get_etat() == ETATZERO) {
return vi ;
}
// Cas general
assert(vi.get_etat() == ETATQCQ) ; // sinon...
if (vi.c == 0x0) { // Il faut la valeur physique
vi.coef_i() ;
}
Valeur vo(vi.get_mg()) ; // Valeur resultat
vo.set_etat_c_qcq() ;
*(vo.c) = abs( *(vi.c) ) ; // abs(Mtbl)
return vo ;
}
//----------------//
// Cube root //
//----------------//
Valeur racine_cubique(const Valeur& vi)
{
// Protection
assert(vi.get_etat() != ETATNONDEF) ;
// Cas ETATZERO
if (vi.get_etat() == ETATZERO) {
return vi ;
}
// Cas general
assert(vi.get_etat() == ETATQCQ) ; // sinon...
if (vi.c == 0x0) { // Il faut la valeur physique
vi.coef_i() ;
}
Valeur vo(vi.get_mg()) ; // Valeur resultat
vo.set_etat_c_qcq() ;
*(vo.c) = racine_cubique( *(vi.c) ) ;
return vo ;
}
//-------------------------------//
// totalmax //
//-------------------------------//
double totalmax(const Valeur& vi) {
// Protection
assert(vi.get_etat() != ETATNONDEF) ;
// Tbl resu(vi.get_mg()->get_nzone()) ;
double resu ;
if (vi.get_etat() == ETATZERO) {
resu = 0 ;
}
else {
assert(vi.get_etat() == ETATQCQ) ;
if (vi.c == 0x0) { // Il faut la valeur physique
vi.coef_i() ;
}
resu = totalmax( *(vi.c) ) ; // max(Mtbl)
}
return resu ;
}
//-------------------------------//
// totalmin //
//-------------------------------//
double totalmin(const Valeur& vi) {
// Protection
assert(vi.get_etat() != ETATNONDEF) ;
double resu ;
if (vi.get_etat() == ETATZERO) {
resu = 0 ;
}
else {
assert(vi.get_etat() == ETATQCQ) ;
if (vi.c == 0x0) { // Il faut la valeur physique
vi.coef_i() ;
}
resu = totalmin( *(vi.c) ) ; // min(Mtbl)
}
return resu ;
}
//-------------------------------//
// max //
//-------------------------------//
Tbl max(const Valeur& vi) {
// Protection
assert(vi.get_etat() != ETATNONDEF) ;
Tbl resu(vi.get_mg()->get_nzone()) ;
if (vi.get_etat() == ETATZERO) {
resu.annule_hard() ;
}
else {
assert(vi.get_etat() == ETATQCQ) ;
if (vi.c == 0x0) { // Il faut la valeur physique
vi.coef_i() ;
}
resu = max( *(vi.c) ) ; // max(Mtbl)
}
return resu ;
}
//-------------------------------//
// min //
//-------------------------------//
Tbl min(const Valeur& vi) {
// Protection
assert(vi.get_etat() != ETATNONDEF) ;
Tbl resu(vi.get_mg()->get_nzone()) ;
if (vi.get_etat() == ETATZERO) {
resu.annule_hard() ;
}
else {
assert(vi.get_etat() == ETATQCQ) ;
if (vi.c == 0x0) { // Il faut la valeur physique
vi.coef_i() ;
}
resu = min( *(vi.c) ) ; // min(Mtbl)
}
return resu ;
}
//-------------------------------//
// norme //
//-------------------------------//
Tbl norme(const Valeur& vi) {
// Protection
assert(vi.get_etat() != ETATNONDEF) ;
Tbl resu(vi.get_mg()->get_nzone()) ;
if (vi.get_etat() == ETATZERO) {
resu.annule_hard() ;
}
else {
assert(vi.get_etat() == ETATQCQ) ;
if (vi.c == 0x0) { // Il faut la valeur physique
vi.coef_i() ;
}
resu = norme( *(vi.c) ) ; // norme(Mtbl)
}
return resu ;
}
//-------------------------------//
// diffrel //
//-------------------------------//
Tbl diffrel(const Valeur& v1, const Valeur& v2) {
// Protections
assert(v1.get_etat() != ETATNONDEF) ;
assert(v2.get_etat() != ETATNONDEF) ;
int nz = v1.get_mg()->get_nzone() ;
Tbl resu(nz) ;
Valeur diff = v1 - v2 ;
Tbl normdiff = norme(diff) ;
Tbl norme2 = norme(v2) ;
assert(normdiff.get_etat() == ETATQCQ) ;
assert(norme2.get_etat() == ETATQCQ) ;
resu.set_etat_qcq() ;
for (int l=0; l<nz; l++) {
if ( norme2(l) == double(0) ) {
resu.set(l) = normdiff(l) ;
}
else{
resu.set(l) = normdiff(l) / norme2(l) ;
}
}
return resu ;
}
//-------------------------------//
// diffrelmax //
//-------------------------------//
Tbl diffrelmax(const Valeur& v1, const Valeur& v2) {
// Protections
assert(v1.get_etat() != ETATNONDEF) ;
assert(v2.get_etat() != ETATNONDEF) ;
int nz = v1.get_mg()->get_nzone() ;
Tbl resu(nz) ;
Tbl max2 = max(abs(v2)) ;
Valeur diff = v1 - v2 ;
Tbl maxdiff = max(abs(diff)) ;
assert(maxdiff.get_etat() == ETATQCQ) ;
assert(max2.get_etat() == ETATQCQ) ;
resu.set_etat_qcq() ;
for (int l=0; l<nz; l++) {
if ( max2(l) == double(0) ) {
resu.set(l) = maxdiff(l) ;
}
else{
resu.set(l) = maxdiff(l) / max2(l) ;
}
}
return resu ;
}
}
Lorene::Mg3d::get_nzone
int get_nzone() const
Returns the number of domains.
Definition grilles.h:448
Lorene::Tbl
Basic array class.
Definition tbl.h:161
Lorene::Tbl::get_etat
int get_etat() const
Gives the logical state.
Definition tbl.h:394
Lorene::Tbl::annule_hard
void annule_hard()
Sets the Tbl to zero in a hard way.
Definition tbl.C:372
Lorene::Tbl::set_etat_qcq
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition tbl.C:361
Lorene::Tbl::set
double & set(int i)
Read/write of a particular element (index i) (1D case).
Definition tbl.h:281
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::Valeur::get_etat
int get_etat() const
Returns the logical state.
Definition valeur.h:726
Lorene::Valeur::c
Mtbl * c
Values of the function at the points of the multi-grid.
Definition valeur.h:299
Lorene::Valeur::coef_i
void coef_i() const
Computes the physical value of *this.
Definition valeur_coef_i.C:140
Lorene::Valeur::get_mg
const Mg3d * get_mg() const
Returns the Mg3d on which the this is defined.
Definition valeur.h:729
Lorene::atan
Cmp atan(const Cmp &)
Arctangent.
Definition cmp_math.C:195
Lorene::sqrt
Cmp sqrt(const Cmp &)
Square root.
Definition cmp_math.C:220
Lorene::log10
Cmp log10(const Cmp &)
Basis 10 logarithm.
Definition cmp_math.C:322
Lorene::exp
Cmp exp(const Cmp &)
Exponential.
Definition cmp_math.C:270
Lorene::sin
Cmp sin(const Cmp &)
Sine.
Definition cmp_math.C:69
Lorene::diffrel
Tbl diffrel(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (norme version).
Definition cmp_math.C:504
Lorene::norme
Tbl norme(const Cmp &)
Sums of the absolute values of all the values of the Cmp in each domain.
Definition cmp_math.C:481
Lorene::acos
Cmp acos(const Cmp &)
Arccosine.
Definition cmp_math.C:169
Lorene::asin
Cmp asin(const Cmp &)
Arcsine.
Definition cmp_math.C:144
Lorene::racine_cubique
Cmp racine_cubique(const Cmp &)
Cube root.
Definition cmp_math.C:245
Lorene::min
Tbl min(const Cmp &)
Minimum values of a Cmp in each domain.
Definition cmp_math.C:458
Lorene::max
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
Definition cmp_math.C:435
Lorene::pow
Cmp pow(const Cmp &, int)
Power .
Definition cmp_math.C:348
Lorene::cos
Cmp cos(const Cmp &)
Cosine.
Definition cmp_math.C:94
Lorene::abs
Cmp abs(const Cmp &)
Absolute value.
Definition cmp_math.C:410
Lorene::diffrelmax
Tbl diffrelmax(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (max version).
Definition cmp_math.C:539
Lorene::tan
Cmp tan(const Cmp &)
Tangent.
Definition cmp_math.C:120
Lorene::log
Cmp log(const Cmp &)
Neperian logarithm.
Definition cmp_math.C:296
Lorene::Heaviside
Mtbl Heaviside(const Mtbl &)
Heaviside function.
Definition mtbl_math.C:317
Lorene::totalmin
double totalmin(const Mtbl &)
Minimum value of the Mtbl elements in all domain.
Definition mtbl_math.C:522
Lorene::totalmax
double totalmax(const Mtbl &)
Maximum value of the Mtbl elements in all domains.
Definition mtbl_math.C:494
Lorene
Lorene prototypes.
Definition app_hor.h:64
Lorene::x
Coord x
x coordinate centered on the grid
Definition map.h:726