31char tensor_calculus_ext_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Tensor/tensor_calculus_ext.C,v 1.14 2014/10/13 08:53:44 j_novak Exp $" ;
109 assert (t1.
mp == t2.
mp) ;
113 Itbl tipe (val_res) ;
115 for (
int i=0 ; i<t1.
valence ; i++)
117 for (
int i=0 ; i<t2.
valence ; i++)
133 Tensor res(*t1.
mp, val_res, tipe, triad_res) ;
138 for (
int i=0 ; i<res.
n_comp ; i++) {
140 for (
int j=0 ; j<t1.
valence ; j++)
141 jeux_indice_t1.
set(j) = jeux_indice_res(j) ;
142 for (
int j=0 ; j<t2.
valence ; j++)
143 jeux_indice_t2.
set(j) = jeux_indice_res(j+t1.
valence) ;
145 res.
set(jeux_indice_res) = t1(jeux_indice_t1)*t2(jeux_indice_t2) ;
167 assert((ind1>=0) && (ind1<val1)) ;
168 assert((ind2>=0) && (ind2<val2)) ;
172 if ( (val1 != 0) && (val2 != 0) ) {
177 int val_res = val1 + val2 - 2;
181 for (
int i=0 ; i<ind1 ; i++)
183 for (
int i=ind1 ; i<val1-1 ; i++)
185 for (
int i=val1-1 ; i<val1+ind2-1 ; i++)
187 for (
int i = val1+ind2-1 ; i<val_res ; i++)
196 Itbl jeux_indice_t1(val1) ;
197 Itbl jeux_indice_t2(val2) ;
203 for (
int j=0 ; j<ind1 ; j++)
204 jeux_indice_t1.
set(j) = jeux_indice_res(j) ;
206 for (
int j=ind1+1 ; j<val1 ; j++)
207 jeux_indice_t1.
set(j) = jeux_indice_res(j-1) ;
209 for (
int j=0 ; j<ind2 ; j++)
210 jeux_indice_t2.
set(j) = jeux_indice_res(val1+j-1) ;
212 for (
int j=ind2+1 ; j<val2 ; j++)
213 jeux_indice_t2.
set(j) = jeux_indice_res(val1+j-2) ;
215 Scalar& work = res.
set(jeux_indice_res) ;
218 for (
int j=1 ; j<=3 ; j++) {
219 jeux_indice_t1.
set(ind1) = j ;
220 jeux_indice_t2.
set(ind2) = j ;
222 work += t1(jeux_indice_t1) % t2(jeux_indice_t2) ;
225 work += t1(jeux_indice_t1) * t2(jeux_indice_t2) ;
239 const Tensor& t2,
int i2,
int j2,
246 assert( val1 >= 2 ) ;
247 assert( val2 >= 2 ) ;
248 assert( (0<=i1) && (i1<j1) && (j1<val1) ) ;
249 assert( (0<=i2) && (i2<j2) && (j2<val2) ) ;
257 int val_res = val1 + val2 - 4 ;
261 for (
int i=0 ; i<i1 ; i++)
264 for (
int i=i1 ; i<j1-1 ; i++)
267 for (
int i=j1-1 ; i<val1-2 ; i++)
270 for (
int i=val1-2 ; i<val1-2+i2 ; i++)
273 for (
int i=val1-2+i2 ; i<val1+j2-3 ; i++)
276 for (
int i=val1+j2-3 ; i<val_res ; i++)
285 Itbl jeux_indice_t1(val1) ;
286 Itbl jeux_indice_t2(val2) ;
292 for (
int k=0 ; k<i1 ; k++)
293 jeux_indice_t1.
set(k) = jeux_indice_res(k) ;
295 for (
int k=i1+1 ; k<j1 ; k++)
296 jeux_indice_t1.
set(k) = jeux_indice_res(k-1) ;
298 for (
int k=j1+1 ; k<val1 ; k++)
299 jeux_indice_t1.
set(k) = jeux_indice_res(k-2) ;
301 for (
int k=0 ; k<i2 ; k++)
302 jeux_indice_t2.
set(k) = jeux_indice_res(val1+k-2) ;
304 for (
int k=i2+1 ; k<j2 ; k++)
305 jeux_indice_t2.
set(k) = jeux_indice_res(val1+k-3) ;
307 for (
int k=j2+1 ; k<val2 ; k++)
308 jeux_indice_t2.
set(k) = jeux_indice_res(val1+k-4) ;
310 Scalar& work = res.
set(jeux_indice_res) ;
313 for (
int i=1 ; i<=3 ; i++) {
315 jeux_indice_t1.
set(i1) = i ;
316 jeux_indice_t2.
set(i2) = i ;
318 for (
int j=1 ; j<=3 ; j++) {
320 jeux_indice_t1.
set(j1) = j ;
321 jeux_indice_t2.
set(j2) = j ;
324 work += t1(jeux_indice_t1) % t2(jeux_indice_t2) ;
327 work += t1(jeux_indice_t1) * t2(jeux_indice_t2) ;
345 assert ((ind_1 >= 0) && (ind_1 < val)) ;
346 assert ((ind_2 >= 0) && (ind_2 < val)) ;
347 assert (ind_1 != ind_2) ;
358 int val_res = val - 2 ;
362 for (
int i=0 ; i<ind_1 ; i++)
364 for (
int i=ind_1 ; i<ind_2-1 ; i++)
366 for (
int i = ind_2-1 ; i<val_res ; i++)
373 Itbl jeux_indice_source(val) ;
379 for (
int j=0 ; j<ind_1 ; j++)
380 jeux_indice_source.
set(j) = jeux_indice_res(j) ;
381 for (
int j=ind_1+1 ; j<ind_2 ; j++)
382 jeux_indice_source.
set(j) = jeux_indice_res(j-1) ;
383 for (
int j=ind_2+1 ; j<val ; j++)
384 jeux_indice_source.
set(j) = jeux_indice_res(j-2) ;
386 Scalar& work = res.
set(jeux_indice_res) ;
389 for (
int j=1 ; j<=3 ; j++) {
390 jeux_indice_source.
set(ind_1) = j ;
391 jeux_indice_source.
set(ind_2) = j ;
392 work += source(jeux_indice_source) ;
411 if (comment != 0x0) ost << comment <<
" : " << endl ;
422 if (n_comp_a >= n_comp_b) {
423 n_comp_max = n_comp_a ;
427 n_comp_max = n_comp_b ;
432 Tbl resu(n_comp_max, nz) ;
437 for (
int ic=0; ic<n_comp_max; ic++) {
441 if (n_comp_max > 1) ost <<
" Comp." ;
442 for (
int j=0 ; j<val ; j++) {
443 ost <<
" " << idx(j) ;
445 if (n_comp_max > 1) ost <<
" : " ;
447 for (
int l=0; l<nz; l++) {
448 ost <<
" " << diff(l) ;
449 resu.
set(ic, l) = diff(l) ;
467 if (comment != 0x0) ost << comment <<
" : " << endl ;
478 if (n_comp_a >= n_comp_b) {
479 n_comp_max = n_comp_a ;
483 n_comp_max = n_comp_b ;
488 Tbl resu(n_comp_max, nz) ;
493 for (
int ic=0; ic<n_comp_max; ic++) {
497 if (n_comp_max > 1) ost <<
" Comp." ;
498 for (
int j=0 ; j<val ; j++) {
499 ost <<
" " << idx(j) ;
501 if (n_comp_max > 1) ost <<
" : " ;
503 for (
int l=0; l<nz; l++) {
504 ost <<
" " << diff(l) ;
505 resu.
set(ic, l) = diff(l) ;
523 if (comment != 0x0) ost << comment <<
" : " << endl ;
530 Tbl resu(n_comp, nz) ;
535 for (
int ic=0; ic<n_comp; ic++) {
538 Tbl diff =
max( aa(idx) ) ;
540 if (val > 0) ost <<
" Comp." ;
541 for (
int j=0 ; j<val ; j++) {
542 ost <<
" " << idx(j) ;
544 if (val > 0) ost <<
" : " ;
546 for (
int l=0; l<nz; l++) {
547 ost <<
" " << diff(l) ;
548 resu.
set(ic, l) = diff(l) ;
566 if (comment != 0x0) ost << comment <<
" : " << endl ;
573 Tbl resu(n_comp, nz) ;
578 for (
int ic=0; ic<n_comp; ic++) {
581 Tbl diff =
min( aa(idx) ) ;
583 if (val > 0) ost <<
" Comp." ;
584 for (
int j=0 ; j<val ; j++) {
585 ost <<
" " << idx(j) ;
587 if (val > 0) ost <<
" : " ;
589 for (
int l=0; l<nz; l++) {
590 ost <<
" " << diff(l) ;
591 resu.
set(ic, l) = diff(l) ;
608 if (comment != 0x0) ost << comment <<
" : " << endl ;
615 Tbl resu(n_comp, nz) ;
620 for (
int ic=0; ic<n_comp; ic++) {
626 if (val > 0) ost <<
" Comp." ;
627 for (
int j=0 ; j<val ; j++) {
628 ost <<
" " << idx(j) ;
630 if (val > 0 ) ost <<
" : " ;
633 for (
int l=0; l<nz; l++) {
634 if (verb) ost <<
" " << diff(l) ;
635 resu.
set(ic, l) = diff(l) ;
637 if (verb) ost <<
"\n" ;
651 if (comment != 0x0) ost << comment <<
" : " << endl ;
661 for (
int ic=0; ic<n_comp; ic++) {
664 double aa_c = aa(idx).val_grid_point(0,0,0,0) ;
665 resu.
set(ic) = aa_c ;
667 if ( comment != 0x0 ) {
668 if ( val > 0 ) ost <<
" Comp." ;
669 for (
int j=0 ; j<val ; j++) {
670 ost <<
" " << idx(j) ;
672 if (val > 0 ) ost <<
" : " ;
674 ost << aa_c << endl ;
690 if (comment != 0x0) ost << comment <<
" : " << endl ;
703 for (
int ic=0; ic<n_comp; ic++) {
706 if (l_excluded != 0) x0 = max_dom(ic, 0) ;
707 else x0 = max_dom(ic, 1) ;
708 for (
int l=0; l<nz; l++) {
709 if (l == l_excluded) continue ;
710 double x = max_dom(ic,l) ;
716 if ( comment != 0x0 ) {
717 if ( val > 0 ) ost <<
" Comp." ;
719 for (
int j=0 ; j<val ; j++) {
720 ost <<
" " << idx(j) ;
722 if (val > 0 ) ost <<
" : " ;
740 if (comment != 0x0) ost << comment <<
" : " << endl ;
753 for (
int ic=0; ic<n_comp; ic++) {
756 if (l_excluded != 0) x0 = min_dom(ic, 0) ;
757 else x0 = min_dom(ic, 1) ;
758 for (
int l=0; l<nz; l++) {
759 if (l == l_excluded) continue ;
760 double x = min_dom(ic,l) ;
766 if ( comment != 0x0 ) {
767 if ( val > 0 ) ost <<
" Comp." ;
769 for (
int j=0 ; j<val ; j++) {
770 ost <<
" " << idx(j) ;
772 if (val > 0 ) ost <<
" : " ;
789 ostream& ost,
bool verb) {
791 if (comment != 0x0) ost << comment <<
" : " << endl ;
793 Tbl maxabs_dom =
maxabs(aa, 0x0, ost, verb) ;
804 for (
int ic=0; ic<n_comp; ic++) {
807 if (l_excluded != 0) x0 = maxabs_dom(ic, 0) ;
808 else x0 = maxabs_dom(ic, 1) ;
809 for (
int l=0; l<nz; l++) {
810 if (l == l_excluded) continue ;
811 double x = maxabs_dom(ic,l) ;
817 if ( comment != 0x0 ) {
818 if ( val > 0 ) ost <<
" Comp." ;
820 for (
int j=0 ; j<val ; j++) {
821 ost <<
" " << idx(j) ;
823 if (val > 0 ) ost <<
" : " ;
Vectorial bases (triads) with respect to which the tensorial components are defined.
Basic integer array class.
int & set(int i)
Read/write of a particular element (index i ) (1D case)
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
int get_nzone() const
Returns the number of domains.
Tensor field of valence 0 (or component of a tensorial field).
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
double & set(int i)
Read/write of a particular element (index i) (1D case)
Base_val operator*(const Base_val &, const Base_val &)
This operator is used when calling multiplication or division of Valeur .
Tbl diffrel(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (norme version).
Tbl min(const Cmp &)
Minimum values of a Cmp in each domain.
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
Cmp abs(const Cmp &)
Absolute value.
Tbl diffrelmax(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (max version).
const Map & get_mp() const
Returns the mapping.
int get_index_type(int i) const
Gives the type (covariant or contravariant) of the index number i .
const Map *const mp
Mapping on which the numerical values at the grid points are defined.
int get_valence() const
Returns the valence.
virtual Itbl indices(int pos) const
Returns the indices of a component given by its position in the array cmp .
int valence
Valence of the tensor (0 = scalar, 1 = vector, etc...)
int get_n_comp() const
Returns the number of stored components.
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
int n_comp
Number of stored components, depending on the symmetry.
Itbl type_indice
1D array of integers (class Itbl ) of size valence containing the type of each index: COV for a cov...
Scalar & set(const Itbl &ind)
Returns the value of a component (read/write version).
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .
Tbl maxabs_all_domains(const Tensor &aa, int l_excluded=-1, const char *comment=0x0, ostream &ost=cout, bool verb=true)
Maximum of the absolute value of each component of a tensor over all the domains.
Tbl maxabs(const Tensor &aa, const char *comment=0x0, ostream &ost=cout, bool verb=true)
Maxima in each domain of the absolute values of the tensor components.
Tbl max_all_domains(const Tensor &aa, int l_excluded=-1, const char *comment=0x0, ostream &ost=cout)
Maximum value of each component of a tensor over all the domains.
Tbl central_value(const Tensor &aa, const char *comment=0x0, ostream &ost=cout)
Central value of each component of a tensor.
Tbl min_all_domains(const Tensor &aa, int l_excluded=-1, const char *comment=0x0, ostream &ost=cout)
Minimum value of each component of a tensor over all the domains.