KADATH
metric_dirac.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 "space.hpp"
21 #include "tensor.hpp"
22 #include "metric.hpp"
23 #include "term_eq.hpp"
24 #include "scalar.hpp"
25 #include "tensor_impl.hpp"
26 #include "system_of_eqs.hpp"
27 #include "metric_tensor.hpp"
28 namespace Kadath {
29 Metric_dirac::Metric_dirac (Metric_tensor& met) :
30  Metric_conf(met) {
31  type_tensor = met.get_type() ;
32 }
33 
34 Metric_dirac::Metric_dirac (const Metric_dirac& so) :
35  Metric_conf (so) {
36 }
37 
38 Metric_dirac::~Metric_dirac() {
39 }
40 
41 
42 
44 
45  // Need christoffels
46  if (p_christo[dd]==0x0)
47  compute_christo(dd) ;
48 
49  // Get the conformal factor
50 /* Tensor res_val (espace, 2, CON, basis) ;
51  Tensor res_der (espace, 2, CON, basis) ;
52 
53 
54  Term_eq der_cov (fmet.derive (COV, ' ', (*p_met_cov[dd]))) ;
55  Term_eq der_con (fmet.derive (COV, ' ', (*p_met_con[dd]))) ;
56  Term_eq dder_con (fmet.derive (COV, ' ', der_con)) ;
57  bool doder = ((der_cov.der_t==0x0) || (der_con.der_t==0x0) || (dder_con.der_t==0x0)) ? false : true ;
58 
59  Index pos (res_val) ;
60  do {
61  Val_domain cmpval (espace.get_domain(dd)) ;
62  cmpval = 0 ;
63 
64  for (int k=1 ; k<=espace.get_ndim() ; k++)
65  for (int l=1 ; l<=espace.get_ndim() ; l++) {
66 
67  cmpval += 0.5 * ((*p_met_con[dd]->val_t)(k,l)(dd)*(*dder_con.val_t)(k,l, pos(0)+1, pos(1)+1)(dd)
68  - (*der_con.val_t)(l, pos(0)+1, k)(dd) * (*der_con.val_t)(k, pos(1)+1, l)(dd)) ;
69  for (int m=1 ; m<=espace.get_ndim() ; m++)
70  for (int n=1 ; n<espace.get_ndim() ; n++)
71  cmpval += 0.5 * (
72  -(*p_met_cov[dd]->val_t)(k,l)(dd)*(*p_met_con[dd]->val_t)(m,n)(dd)*(*der_con.val_t)(m, pos(0)+1, k)(dd)*(*der_con.val_t)(n, pos(1)+1, l)(dd)
73  +(*p_met_cov[dd]->val_t)(m,l)(dd)*(*p_met_con[dd]->val_t)(pos(0)+1,k)(dd)*(*der_con.val_t)(k, m, n)(dd)*(*der_con.val_t)(n, pos(1)+1, l)(dd)
74  +(*p_met_cov[dd]->val_t)(k,n)(dd)*(*p_met_con[dd]->val_t)(pos(1)+1,l)(dd)*(*der_con.val_t)(l, m, n)(dd)*(*der_con.val_t)(m, pos(0)+1, k)(dd)
75  +0.5*(*p_met_cov[dd]->val_t)(pos(0)+1,k)(dd)*(*p_met_con[dd]->val_t)(pos(1)+1,l)(dd)*(*der_cov.val_t)(k, m, n)(dd)*(*der_con.val_t)(l,m,n)(dd)) ;
76 
77  }
78 
79  res_val.set(pos).set_domain(dd) = cmpval ;
80 
81  if (doder) {
82  Val_domain cmpder (espace.get_domain(dd)) ;
83  cmpder = 0 ;
84 
85  for (int k=1 ; k<=espace.get_ndim() ; k++)
86  for (int l=1 ; l<=espace.get_ndim() ; l++) {
87 
88  cmpder += 0.5 * ((*p_met_con[dd]->der_t)(k,l)(dd)*(*dder_con.val_t)(k,l, pos(0)+1, pos(1)+1)(dd)
89  +(*p_met_con[dd]->val_t)(k,l)(dd)*(*dder_con.der_t)(k,l, pos(0)+1, pos(1)+1)(dd)
90  - (*der_con.der_t)(l, pos(0)+1, k)(dd) * (*der_con.val_t)(k, pos(1)+1, l)(dd)
91  - (*der_con.val_t)(l, pos(0)+1, k)(dd) * (*der_con.der_t)(k, pos(1)+1, l)(dd)) ;
92  for (int m=1 ; m<=espace.get_ndim() ; m++)
93  for (int n=1 ; n<espace.get_ndim() ; n++)
94  cmpder += 0.5 * (
95  -(*p_met_cov[dd]->der_t)(k,l)(dd)*(*p_met_con[dd]->val_t)(m,n)(dd)*(*der_con.val_t)(m, pos(0)+1, k)(dd)*(*der_con.val_t)(n, pos(1)+1, l)(dd)
96  -(*p_met_cov[dd]->val_t)(k,l)(dd)*(*p_met_con[dd]->der_t)(m,n)(dd)*(*der_con.val_t)(m, pos(0)+1, k)(dd)*(*der_con.val_t)(n, pos(1)+1, l)(dd)
97  -(*p_met_cov[dd]->val_t)(k,l)(dd)*(*p_met_con[dd]->val_t)(m,n)(dd)*(*der_con.der_t)(m, pos(0)+1, k)(dd)*(*der_con.val_t)(n, pos(1)+1, l)(dd)
98  -(*p_met_cov[dd]->val_t)(k,l)(dd)*(*p_met_con[dd]->val_t)(m,n)(dd)*(*der_con.val_t)(m, pos(0)+1, k)(dd)*(*der_con.der_t)(n, pos(1)+1, l)(dd)
99  +(*p_met_cov[dd]->der_t)(m,l)(dd)*(*p_met_con[dd]->val_t)(pos(0)+1,k)(dd)*(*der_con.val_t)(k, m, n)(dd)*(*der_con.val_t)(n, pos(1)+1, l)(dd)
100  +(*p_met_cov[dd]->val_t)(m,l)(dd)*(*p_met_con[dd]->der_t)(pos(0)+1,k)(dd)*(*der_con.val_t)(k, m, n)(dd)*(*der_con.val_t)(n, pos(1)+1, l)(dd)
101  +(*p_met_cov[dd]->val_t)(m,l)(dd)*(*p_met_con[dd]->val_t)(pos(0)+1,k)(dd)*(*der_con.der_t)(k, m, n)(dd)*(*der_con.val_t)(n, pos(1)+1, l)(dd)
102  +(*p_met_cov[dd]->val_t)(m,l)(dd)*(*p_met_con[dd]->val_t)(pos(0)+1,k)(dd)*(*der_con.val_t)(k, m, n)(dd)*(*der_con.der_t)(n, pos(1)+1, l)(dd)
103  +(*p_met_cov[dd]->der_t)(k,n)(dd)*(*p_met_con[dd]->val_t)(pos(1)+1,l)(dd)*(*der_con.val_t)(l, m, n)(dd)*(*der_con.val_t)(m, pos(0)+1, k)(dd)
104  +(*p_met_cov[dd]->val_t)(k,n)(dd)*(*p_met_con[dd]->der_t)(pos(1)+1,l)(dd)*(*der_con.val_t)(l, m, n)(dd)*(*der_con.val_t)(m, pos(0)+1, k)(dd)
105  +(*p_met_cov[dd]->val_t)(k,n)(dd)*(*p_met_con[dd]->val_t)(pos(1)+1,l)(dd)*(*der_con.der_t)(l, m, n)(dd)*(*der_con.val_t)(m, pos(0)+1, k)(dd)
106  +(*p_met_cov[dd]->val_t)(k,n)(dd)*(*p_met_con[dd]->val_t)(pos(1)+1,l)(dd)*(*der_con.val_t)(l, m, n)(dd)*(*der_con.der_t)(m, pos(0)+1, k)(dd)
107  +0.5*(*p_met_cov[dd]->der_t)(pos(0)+1,k)(dd)*(*p_met_con[dd]->val_t)(pos(1)+1,l)(dd)*(*der_cov.val_t)(k, m, n)(dd)*(*der_con.val_t)(l,m,n)(dd)
108  +0.5*(*p_met_cov[dd]->val_t)(pos(0)+1,k)(dd)*(*p_met_con[dd]->der_t)(pos(1)+1,l)(dd)*(*der_cov.val_t)(k, m, n)(dd)*(*der_con.val_t)(l,m,n)(dd)
109  +0.5*(*p_met_cov[dd]->val_t)(pos(0)+1,k)(dd)*(*p_met_con[dd]->val_t)(pos(1)+1,l)(dd)*(*der_cov.der_t)(k, m, n)(dd)*(*der_con.val_t)(l,m,n)(dd)
110  +0.5*(*p_met_cov[dd]->val_t)(pos(0)+1,k)(dd)*(*p_met_con[dd]->val_t)(pos(1)+1,l)(dd)*(*der_cov.val_t)(k, m, n)(dd)*(*der_con.der_t)(l,m,n)(dd)
111  ) ;
112 
113 
114  }
115 
116  res_der.set(pos).set_domain(dd) = cmpder ;
117  }
118  }
119 
120  while (pos.inc()) ;
121 
122 */
123 
124  Tensor res_val (espace, 2, COV, basis) ;
125  Tensor res_der (espace, 2, COV, basis) ;
126 
127 
128  Term_eq der_cov (fmet.derive (COV, ' ', (*p_met_cov[dd]))) ;
129  Term_eq dder_cov (fmet.derive (COV, ' ', der_cov)) ;
130  Term_eq der_con (fmet.derive (COV, ' ', (*p_met_con[dd]))) ;
131  bool doder = ((der_cov.der_t==0x0) || (der_con.der_t==0x0) || (dder_cov.der_t==0x0)) ? false : true ;
132 
133  Index pos (res_val) ;
134  do {
135  Val_domain cmpval (espace.get_domain(dd)) ;
136  cmpval = 0 ;
137 
138  for (int k=1 ; k<=espace.get_ndim() ; k++)
139  for (int l=1 ; l<=espace.get_ndim() ; l++) {
140 
141  cmpval += - 0.5 * ((*p_met_con[dd]->val_t)(k,l)(dd)*(*dder_cov.val_t)(k,l, pos(0)+1, pos(1)+1)(dd)
142  +(*der_con.val_t)(pos(0)+1, k, l)(dd) * (*der_cov.val_t)(k, pos(1)+1, l)(dd)
143  + (*der_con.val_t)(pos(1)+1, k, l)(dd) * (*der_cov.val_t)(k, pos(0)+1, l)(dd))
144  - (*p_christo[dd]->val_t)(pos(0)+1,l,k)(dd)*(*p_christo[dd]->val_t)(pos(1)+1,k,l)(dd) ;
145  }
146 
147  res_val.set(pos).set_domain(dd) = cmpval ;
148 
149  if (doder) {
150  Val_domain cmpder (espace.get_domain(dd)) ;
151  cmpder = 0 ;
152 
153  for (int k=1 ; k<=espace.get_ndim() ; k++)
154  for (int l=1 ; l<=espace.get_ndim() ; l++) {
155 
156  cmpder += - 0.5 * ((*p_met_con[dd]->der_t)(k,l)(dd)*(*dder_cov.val_t)(k,l, pos(0)+1, pos(1)+1)(dd)
157  + (*p_met_con[dd]->val_t)(k,l)(dd)*(*dder_cov.der_t)(k,l, pos(0)+1, pos(1)+1)(dd)
158  +(*der_con.der_t)(pos(0)+1, k, l)(dd) * (*der_cov.val_t)(k, pos(1)+1, l)(dd)
159  +(*der_con.val_t)(pos(0)+1, k, l)(dd) * (*der_cov.der_t)(k, pos(1)+1, l)(dd)
160  + (*der_con.der_t)(pos(1)+1, k, l)(dd) * (*der_cov.val_t)(k, pos(0)+1, l)(dd)
161  + (*der_con.val_t)(pos(1)+1, k, l)(dd) * (*der_cov.der_t)(k, pos(0)+1, l)(dd))
162  - (*p_christo[dd]->der_t)(pos(0)+1,l,k)(dd)*(*p_christo[dd]->val_t)(pos(1)+1,k,l)(dd)
163  - (*p_christo[dd]->val_t)(pos(0)+1,l,k)(dd)*(*p_christo[dd]->der_t)(pos(1)+1,k,l)(dd) ;
164  }
165 
166  res_der.set(pos).set_domain(dd) = cmpder ;
167  }
168  }
169 
170  while (pos.inc()) ;
171 
172  if (!doder) {
173  if (p_ricci_tensor[dd]==0x0)
174  p_ricci_tensor[dd] = new Term_eq(dd, res_val) ;
175  else
176  *p_ricci_tensor[dd] = Term_eq (dd, res_val) ;
177  }
178  else {
179  if (p_ricci_tensor[dd]==0x0)
180  p_ricci_tensor[dd] = new Term_eq(dd, res_val, res_der) ;
181  else
182  *p_ricci_tensor[dd] = Term_eq(dd, res_val, res_der) ;
183  }
184 
185 
186 }
187 
188 
190 
191  // Need that
192  if (p_met_con[dd]==0x0)
193  compute_con(dd) ;
194  if (p_ricci_tensor[dd]==0x0)
196 
197  bool doder = ((p_met_con[dd]->der_t==0x0) || (p_ricci_tensor[dd]->der_t==0x0)) ? false : true ;
198  Scalar res_val (espace) ;
199  Scalar res_der (espace) ;
200 
201  Val_domain cmpval (espace.get_domain(dd)) ;
202  cmpval = 0 ;
203 
204  for (int i=1 ; i<=espace.get_ndim() ; i++)
205  for (int j=1 ; j<=espace.get_ndim() ; j++)
206  cmpval += (*p_met_con[dd]->val_t)(i,j)(dd)*(*p_ricci_tensor[dd]->val_t)(i,j)(dd) ;
207  res_val.set_domain(dd) = cmpval ;
208 
209  if (doder) {
210  Val_domain cmpder (espace.get_domain(dd)) ;
211  cmpder = 0 ;
212 
213  for (int i=1 ; i<=espace.get_ndim() ; i++)
214  for (int j=1 ; j<=espace.get_ndim() ; j++)
215  cmpder += (*p_met_con[dd]->val_t)(i,j)(dd)*(*p_ricci_tensor[dd]->der_t)(i,j)(dd)
216  + (*p_met_con[dd]->der_t)(i,j)(dd)*(*p_ricci_tensor[dd]->val_t)(i,j)(dd) ;
217  res_der.set_domain(dd) = cmpder ;
218  }
219 
220  if (!doder) {
221  if (p_ricci_scalar[dd]==0x0)
222  p_ricci_scalar[dd] = new Term_eq(dd, res_val) ;
223  else
224  *p_ricci_scalar[dd] = Term_eq (dd, res_val) ;
225  }
226  else {
227  if (p_ricci_scalar[dd]==0x0)
228  p_ricci_scalar[dd] = new Term_eq(dd, res_val, res_der) ;
229  else
230  *p_ricci_scalar[dd] = Term_eq(dd, res_val, res_der) ;
231  }
232 }}
Class that gives the position inside a multi-dimensional Array.
Definition: index.hpp:38
bool inc(int increm, int var=0)
Increments the position of the Index.
Definition: index.hpp:99
Class to deal with a metric which determinant is 1.
Definition: metric.hpp:443
virtual void compute_con(int) const
Computes the contravariant representation, in a given Domain.
Metric_flat fmet
Associated flat metric.
Definition: metric.hpp:448
const Base_tensor & basis
The tensorial basis used.
Definition: metric.hpp:447
virtual void compute_christo(int) const
Computes the Christoffel symbols, in a given Domain.
Class to deal with a conformal metric in the Dirac gauge.
Definition: metric.hpp:485
virtual void compute_ricci_tensor(int) const
Computes the Ricci tensor, in a given Domain.
virtual void compute_ricci_scalar(int) const
Computes the Ricci scalar, in a given Domain.
virtual Term_eq derive(int, char, const Term_eq &) const
Computes the covariant derivative of a Term_eq (assumes Cartesian basis of decomposition).
Particular type of Tensor, dedicated to the desription of metrics.
int type_tensor
States if one works in the CON or COV representation.
Definition: metric.hpp:86
MMPtr_array< Term_eq > p_ricci_tensor
Array of pointers on various Term_eq.
Definition: metric.hpp:70
const Space & espace
The associated Space.
Definition: metric.hpp:42
MMPtr_array< Term_eq > p_met_cov
Array of pointers on various Term_eq.
Definition: metric.hpp:50
MMPtr_array< Term_eq > p_met_con
Array of pointers on various Term_eq.
Definition: metric.hpp:55
MMPtr_array< Term_eq > p_christo
Array of pointers on various Term_eq.
Definition: metric.hpp:60
MMPtr_array< Term_eq > p_ricci_scalar
Array of pointers on various Term_eq.
Definition: metric.hpp:75
The class Scalar does not really implements scalars in the mathematical sense but rather tensorial co...
Definition: scalar.hpp:67
Val_domain & set_domain(int)
Read/write of a particular Val_domain.
Definition: scalar.hpp:555
const Domain * get_domain(int i) const
returns a pointer on the domain.
Definition: space.hpp:1385
int get_ndim() const
Returns the number of dimensions.
Definition: space.hpp:1373
Tensor handling.
Definition: tensor.hpp:149
Scalar & set(const Array< int > &ind)
Returns the value of a component (read/write version).
Definition: tensor_impl.hpp:91
This class is intended to describe the manage objects appearing in the equations.
Definition: term_eq.hpp:62
Tensor * der_t
Pointer on the variation, if the Term_eq is a Tensor.
Definition: term_eq.hpp:69
Tensor * val_t
Pointer on the value, if the Term_eq is a Tensor.
Definition: term_eq.hpp:68
Class for storing the basis of decompositions of a field and its values on both the configuration and...
Definition: val_domain.hpp:69