KADATH
ope_determinant.cpp
1 /*
2  Copyright 2017 Philippe Grandclement
3 
4  This file is part of Kadath.
5 
6  Kadath is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 3 of the License, or
9  (at your option) any later version.
10 
11  Kadath is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with Kadath. If not, see <http://www.gnu.org/licenses/>.
18 */
19 
20 #include "ope_eq.hpp"
21 #include "scalar.hpp"
22 #include "tensor_impl.hpp"
23 namespace Kadath {
24 Ope_determinant::Ope_determinant (const System_of_eqs* zesys, Ope_eq* target) : Ope_eq(zesys, target->get_dom(), 1) {
25  parts[0] = target ;
26 }
27 
28 Ope_determinant::~Ope_determinant() {
29 }
30 
31 Term_eq Ope_determinant::action() const {
32 
33  Term_eq target (parts[0]->action()) ;
34 
35  // Various checks
36  if (target.type_data != TERM_T) {
37  cerr << "Ope_determinant only defined with respect for a tensor" << endl ;
38  abort() ;
39  }
40  if (target.val_t->get_valence()!=2) {
41  cerr << "Ope_determinant only defined with respect to a second order tensor" << endl ;
42  abort() ;
43  }
44  if (target.val_t->get_index_type(0)!=target.val_t->get_index_type(1)) {
45  cerr << "Ope_determinant only defined with respect to a tensor which indices are of the same type" << endl ;
46  abort() ;
47  }
48 
49  int ndim = target.val_t->get_space().get_ndim() ;
50 
51  switch (ndim) {
52  case 2 :{
53 
54  Val_domain resval ((*target.val_t)(1,1)(dom)* (*target.val_t)(2,2)(dom) -
55  (*target.val_t)(1,2)(dom)* (*target.val_t)(2,1)(dom)) ;
56  Scalar res (target.val_t->get_space()) ;
57  res.set_domain(dom) = resval ;
58 
59  if (target.der_t != 0x0) {
60  Val_domain derval ((*target.val_t)(1,1)(dom)* (*target.der_t)(2,2)(dom) +
61  (*target.der_t)(1,1)(dom)* (*target.val_t)(2,2)(dom) -
62  (*target.val_t)(1,2)(dom)* (*target.der_t)(2,1)(dom) -
63  (*target.der_t)(1,2)(dom)* (*target.val_t)(2,1)(dom)) ;
64 
65  Scalar der (target.val_t->get_space()) ;
66  der.set_domain(dom) = derval ;
67 
68  return Term_eq (dom, res, der) ;
69  }
70  else {
71  return Term_eq (dom, res) ;
72  }
73  }
74  case 3 : {
75 
76 
77  Val_domain resval ((*target.val_t)(1,1)(dom)* (*target.val_t)(2,2)(dom)* (*target.val_t)(3,3)(dom)
78  + (*target.val_t)(2,1)(dom)* (*target.val_t)(3,2)(dom)* (*target.val_t)(1,3)(dom)
79  + (*target.val_t)(3,1)(dom)* (*target.val_t)(1,2)(dom)* (*target.val_t)(2,3)(dom)
80  - (*target.val_t)(1,3)(dom)* (*target.val_t)(2,2)(dom)* (*target.val_t)(3,1)(dom)
81  - (*target.val_t)(2,3)(dom)* (*target.val_t)(3,2)(dom)* (*target.val_t)(1,1)(dom)
82  - (*target.val_t)(3,3)(dom)* (*target.val_t)(1,2)(dom)* (*target.val_t)(2,1)(dom)) ;
83  Scalar res (target.val_t->get_space()) ;
84  res.set_domain(dom) = resval ;
85 
86  //if (target.val_t->get_basis().get_basis(dom)==SPHERICAL_BASIS) {
87  // res.affect_parameters() ;
88  // res.set_parameters()->set_m_order() = 2 ;
89  //}
90 
91  if (target.der_t!=0x0) {
92  Val_domain derval ((*target.der_t)(1,1)(dom)* (*target.val_t)(2,2)(dom)* (*target.val_t)(3,3)(dom)
93  + (*target.val_t)(1,1)(dom)* (*target.der_t)(2,2)(dom)* (*target.val_t)(3,3)(dom)
94  + (*target.val_t)(1,1)(dom)* (*target.val_t)(2,2)(dom)* (*target.der_t)(3,3)(dom)
95  + (*target.der_t)(2,1)(dom)* (*target.val_t)(3,2)(dom)* (*target.val_t)(1,3)(dom)
96  + (*target.val_t)(2,1)(dom)* (*target.der_t)(3,2)(dom)* (*target.val_t)(1,3)(dom)
97  + (*target.val_t)(2,1)(dom)* (*target.val_t)(3,2)(dom)* (*target.der_t)(1,3)(dom)
98  + (*target.der_t)(3,1)(dom)* (*target.val_t)(1,2)(dom)* (*target.val_t)(2,3)(dom)
99  + (*target.val_t)(3,1)(dom)* (*target.der_t)(1,2)(dom)* (*target.val_t)(2,3)(dom)
100  + (*target.val_t)(3,1)(dom)* (*target.val_t)(1,2)(dom)* (*target.der_t)(2,3)(dom)
101  - (*target.der_t)(1,3)(dom)* (*target.val_t)(2,2)(dom)* (*target.val_t)(3,1)(dom)
102  - (*target.val_t)(1,3)(dom)* (*target.der_t)(2,2)(dom)* (*target.val_t)(3,1)(dom)
103  - (*target.val_t)(1,3)(dom)* (*target.val_t)(2,2)(dom)* (*target.der_t)(3,1)(dom)
104  - (*target.der_t)(2,3)(dom)* (*target.val_t)(3,2)(dom)* (*target.val_t)(1,1)(dom)
105  - (*target.val_t)(2,3)(dom)* (*target.der_t)(3,2)(dom)* (*target.val_t)(1,1)(dom)
106  - (*target.val_t)(2,3)(dom)* (*target.val_t)(3,2)(dom)* (*target.der_t)(1,1)(dom)
107  - (*target.der_t)(3,3)(dom)* (*target.val_t)(1,2)(dom)* (*target.val_t)(2,1)(dom)
108  - (*target.val_t)(3,3)(dom)* (*target.der_t)(1,2)(dom)* (*target.val_t)(2,1)(dom)
109  - (*target.val_t)(3,3)(dom)* (*target.val_t)(1,2)(dom)* (*target.der_t)(2,1)(dom)) ;
110 
111  Scalar der (target.val_t->get_space()) ;
112  der.set_domain(dom) = derval ;
113 
114  // if (target.val_t->get_basis().get_basis(dom)==SPHERICAL_BASIS) {
115  // der.affect_parameters() ;
116  // der.set_parameters()->set_m_order() = 2 ;
118  return Term_eq (dom, res, der) ;
119  }
120  else
121  return Term_eq (dom, res) ;
122  }
123  default:
124  cerr << "Ope_determinant not implemented in " << ndim << " dimensions" << endl ;
125  abort() ;
126  }
127 }}
Kadath::Ope_determinant::~Ope_determinant
~Ope_determinant() override
Destructor.
Definition: ope_determinant.cpp:28
Kadath::Ope_determinant::action
Term_eq action() const override
Computes the action of the current Ope_eq using its various parts.
Definition: ope_determinant.cpp:31
Kadath::Ope_determinant::Ope_determinant
Ope_determinant(const System_of_eqs *syst, Ope_eq *so)
Constructor.
Definition: ope_determinant.cpp:24
Kadath::Ope_eq
Abstract class that describes the various operators that can appear in the equations.
Definition: ope_eq.hpp:32
Kadath::Ope_eq::parts
MMPtr_array< Ope_eq > parts
Pointers of the various parts of the current operator.
Definition: ope_eq.hpp:38
Kadath::Ope_eq::dom
int dom
Index of the Domain where the operator is defined.
Definition: ope_eq.hpp:36
Kadath::Scalar
The class Scalar does not really implements scalars in the mathematical sense but rather tensorial co...
Definition: scalar.hpp:67
Kadath::Scalar::set_domain
Val_domain & set_domain(int)
Read/write of a particular Val_domain.
Definition: scalar.hpp:555
Kadath::Space::get_ndim
int get_ndim() const
Returns the number of dimensions.
Definition: space.hpp:1373
Kadath::System_of_eqs
Class used to describe and solve a system of equations.
Definition: system_of_eqs.hpp:60
Kadath::Tensor::get_index_type
int get_index_type(int i) const
Gives the type (covariant or contravariant) of a given index.
Definition: tensor.hpp:526
Kadath::Tensor::get_valence
int get_valence() const
Returns the valence.
Definition: tensor.hpp:509
Kadath::Tensor::get_space
const Space & get_space() const
Returns the Space.
Definition: tensor.hpp:499
Kadath::Term_eq
This class is intended to describe the manage objects appearing in the equations.
Definition: term_eq.hpp:62
Kadath::Term_eq::der_t
Tensor * der_t
Pointer on the variation, if the Term_eq is a Tensor.
Definition: term_eq.hpp:69
Kadath::Term_eq::type_data
const int type_data
Flag describing the type of data :
Definition: term_eq.hpp:75
Kadath::Term_eq::val_t
Tensor * val_t
Pointer on the value, if the Term_eq is a Tensor.
Definition: term_eq.hpp:68
Kadath::Val_domain
Class for storing the basis of decompositions of a field and its values on both the configuration and...
Definition: val_domain.hpp:69