LORENE
time_slice.C
1 /*
2  * Methods of class Time_slice
3  *
4  * (see file time_slice.h for documentation).
5  *
6  */
7 
8 /*
9  * Copyright (c) 2004 Eric Gourgoulhon, Jose Luis Jaramillo & Jerome Novak
10  *
11  * This file is part of LORENE.
12  *
13  * LORENE is free software; you can redistribute it and/or modify
14  * it under the terms of the GNU General Public License version 2
15  * as published by the Free Software Foundation.
16  *
17  * LORENE is distributed in the hope that it will be useful,
18  * but WITHOUT ANY WARRANTY; without even the implied warranty of
19  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20  * GNU General Public License for more details.
21  *
22  * You should have received a copy of the GNU General Public License
23  * along with LORENE; if not, write to the Free Software
24  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
25  *
26  */
27 
28 char time_slice_C[] = "$Header: /cvsroot/Lorene/C++/Source/Time_slice/time_slice.C,v 1.16 2014/10/13 08:53:47 j_novak Exp $" ;
29 
30 /*
31  * $Id: time_slice.C,v 1.16 2014/10/13 08:53:47 j_novak Exp $
32  * $Log: time_slice.C,v $
33  * Revision 1.16 2014/10/13 08:53:47 j_novak
34  * Lorene classes and functions now belong to the namespace Lorene.
35  *
36  * Revision 1.15 2014/10/06 15:13:21 j_novak
37  * Modified #include directives to use c++ syntax.
38  *
39  * Revision 1.14 2008/12/02 15:02:21 j_novak
40  * Implementation of the new constrained formalism, following Cordero et al. 2009
41  * paper. The evolution eqs. are solved as a first-order system. Not tested yet!
42  *
43  * Revision 1.13 2004/05/31 09:08:18 e_gourgoulhon
44  * Method sauve and constructor from binary file are now operational.
45  *
46  * Revision 1.12 2004/05/27 15:25:04 e_gourgoulhon
47  * Added constructors from binary file, as well as corresponding
48  * functions sauve and save.
49  *
50  * Revision 1.11 2004/05/12 15:24:20 e_gourgoulhon
51  * Reorganized the #include 's, taking into account that
52  * time_slice.h contains now an #include "metric.h".
53  *
54  * Revision 1.10 2004/05/10 09:08:34 e_gourgoulhon
55  * Added "adm_mass_evol.downdate(jtime)" in method del_deriv.
56  * Added printing of ADM mass in operator>>(ostream&).
57  *
58  * Revision 1.9 2004/05/09 20:57:34 e_gourgoulhon
59  * Added data member adm_mass_evol.
60  *
61  * Revision 1.8 2004/05/05 14:26:25 e_gourgoulhon
62  * Minor modif. in operator>>(ostream& ).
63  *
64  * Revision 1.7 2004/04/07 07:58:21 e_gourgoulhon
65  * Constructor as Minkowski slice: added call to std_spectral_base()
66  * after setting the lapse to 1.
67  *
68  * Revision 1.6 2004/04/01 16:09:02 j_novak
69  * Trace of K_ij is now member of Time_slice (it was member of Time_slice_conf).
70  * Added new methods for checking 3+1 Einstein equations (preliminary).
71  *
72  * Revision 1.5 2004/03/29 11:59:23 e_gourgoulhon
73  * Added operator>>.
74  *
75  * Revision 1.4 2004/03/28 21:29:45 e_gourgoulhon
76  * Evolution_std's renamed with suffix "_evol"
77  * Method gam() modified
78  * Added special constructor for derived classes.
79  *
80  * Revision 1.3 2004/03/26 13:33:02 j_novak
81  * New methods for accessing/updating members (nn(), beta(), gam_uu(), k_uu(), ...)
82  *
83  * Revision 1.2 2004/03/26 08:22:56 e_gourgoulhon
84  * Modifications to take into account the new setting of class
85  * Evolution.
86  *
87  * Revision 1.1 2004/03/24 14:57:17 e_gourgoulhon
88  * First version
89  *
90  *
91  * $Header: /cvsroot/Lorene/C++/Source/Time_slice/time_slice.C,v 1.16 2014/10/13 08:53:47 j_novak Exp $
92  *
93  */
94 
95 // C headers
96 #include <cassert>
97 
98 // Lorene headers
99 #include "time_slice.h"
100 #include "utilitaires.h"
101 
102 
103 
104  //--------------//
105  // Constructors //
106  //--------------//
107 
108 
109 namespace Lorene {
110 // Standard constructor (Hamiltonian-like)
111 // ---------------------------------------
112 Time_slice::Time_slice(const Scalar& lapse_in, const Vector& shift_in,
113  const Sym_tensor& gamma_in, const Sym_tensor& kk_in,
114  int depth_in)
115  : depth(depth_in),
116  scheme_order(depth_in-1),
117  jtime(0),
118  the_time(0., depth_in),
119  gam_dd_evol(depth_in),
120  gam_uu_evol(depth_in),
121  k_dd_evol(depth_in),
122  k_uu_evol(depth_in),
123  n_evol(lapse_in, depth_in),
124  beta_evol(shift_in, depth_in),
125  trk_evol(depth_in),
126  adm_mass_evol() {
127 
128  set_der_0x0() ;
129 
130  double time_init = the_time[jtime] ;
131 
132  p_gamma = new Metric(gamma_in) ;
133 
134  if (gamma_in.get_index_type(0) == COV) {
135  gam_dd_evol.update(gamma_in, jtime, time_init) ;
136  }
137  else {
138  gam_uu_evol.update(gamma_in, jtime, time_init) ;
139  }
140 
141  if (kk_in.get_index_type(0) == COV) {
142  k_dd_evol.update(kk_in, jtime, time_init) ;
143  }
144  else {
145  k_uu_evol.update(kk_in, jtime, time_init) ;
146  }
147 
148  trk_evol.update( kk_in.trace(*p_gamma), jtime, the_time[jtime] ) ;
149 
150 }
151 
152 
153 // Standard constructor (Lagrangian-like)
154 // ---------------------------------------
155 Time_slice::Time_slice(const Scalar& lapse_in, const Vector& shift_in,
156  const Evolution_std<Sym_tensor>& gamma_in)
157  : depth(gamma_in.get_size()),
158  scheme_order(gamma_in.get_size()-1),
159  jtime(0),
160  the_time(0., gamma_in.get_size()),
161  gam_dd_evol( gamma_in.get_size() ),
162  gam_uu_evol( gamma_in.get_size() ),
163  k_dd_evol( gamma_in.get_size() ),
164  k_uu_evol( gamma_in.get_size() ),
165  n_evol(lapse_in, gamma_in.get_size() ),
166  beta_evol(shift_in, gamma_in.get_size() ),
167  trk_evol(gamma_in.get_size() ),
168  adm_mass_evol() {
169 
170  cerr <<
171  "Time_slice constuctor from evolution of gamma not implemented yet !\n" ;
172  abort() ;
173 
174  set_der_0x0() ;
175 
176 }
177 
178 // Constructor as standard time slice of flat spacetime (Minkowski)
179 //-----------------------------------------------------------------
180 Time_slice::Time_slice(const Map& mp, const Base_vect& triad, int depth_in)
181  : depth(depth_in),
182  scheme_order(depth_in-1),
183  jtime(0),
184  the_time(0., depth_in),
185  gam_dd_evol(depth_in),
186  gam_uu_evol(depth_in),
187  k_dd_evol(depth_in),
188  k_uu_evol(depth_in),
189  n_evol(depth_in),
190  beta_evol(depth_in),
191  trk_evol(depth_in),
192  adm_mass_evol() {
193 
194  double time_init = the_time[jtime] ;
195 
196  const Base_vect_spher* ptriad_s =
197  dynamic_cast<const Base_vect_spher*>(&triad) ;
198  bool spher = (ptriad_s != 0x0) ;
199 
200  if (spher) {
201  gam_dd_evol.update( mp.flat_met_spher().cov(), jtime, time_init) ;
202  }
203  else {
204  assert( dynamic_cast<const Base_vect_cart*>(&triad) != 0x0) ;
205  gam_dd_evol.update( mp.flat_met_cart().cov(), jtime, time_init) ;
206  }
207 
208 
209  // K_ij identically zero:
210  Sym_tensor ktmp(mp, COV, triad) ;
211  ktmp.set_etat_zero() ;
212  k_dd_evol.update(ktmp, jtime, time_init) ;
213 
214  // Lapse identically one:
215  Scalar tmp(mp) ;
216  tmp.set_etat_one() ;
217  tmp.std_spectral_base() ;
218  n_evol.update(tmp, jtime, time_init) ;
219 
220  // shift identically zero:
221  Vector btmp(mp, CON, triad) ;
222  btmp.set_etat_zero() ;
223  beta_evol.update(btmp, jtime, time_init) ;
224 
225  // trace(K) identically zero:
226  tmp.set_etat_zero() ;
227  trk_evol.update(tmp, jtime, time_init) ;
228 
229  set_der_0x0() ;
230 }
231 
232 // Constructor from binary file
233 // ----------------------------
234 
235 Time_slice::Time_slice(const Map& mp, const Base_vect& triad, FILE* fich,
236  bool partial_read, int depth_in)
237  : depth(depth_in),
238  the_time(depth_in),
239  gam_dd_evol(depth_in),
240  gam_uu_evol(depth_in),
241  k_dd_evol(depth_in),
242  k_uu_evol(depth_in),
243  n_evol(depth_in),
244  beta_evol(depth_in),
245  trk_evol(depth_in),
246  adm_mass_evol() {
247 
248  // Reading various integer parameters
249  // ----------------------------------
250 
251  int depth_file ;
252  fread_be(&depth_file, sizeof(int), 1, fich) ;
253  if (depth_file != depth_in) {
254  cout <<
255  "Time_slice constructor from file: the depth read in file \n"
256  << " is different from that given in the argument list : \n"
257  << " depth_file = " << depth_file
258  << " <-> depth_in " << depth_in << " !" << endl ;
259  abort() ;
260  }
261  fread_be(&scheme_order, sizeof(int), 1, fich) ;
262  fread_be(&jtime, sizeof(int), 1, fich) ;
263 
264  // Reading the_time
265  // ----------------
266  int jmin = jtime - depth + 1 ;
267  int indicator ;
268  for (int j=jmin; j<=jtime; j++) {
269  fread_be(&indicator, sizeof(int), 1, fich) ;
270  if (indicator == 1) {
271  double xx ;
272  fread_be(&xx, sizeof(double), 1, fich) ;
273  the_time.update(xx, j, xx) ;
274  }
275  }
276 
277  // Reading of various fields
278  // -------------------------
279 
280  // N
281  for (int j=jmin; j<=jtime; j++) {
282  fread_be(&indicator, sizeof(int), 1, fich) ;
283  if (indicator == 1) {
284  Scalar nn_file(mp, *(mp.get_mg()), fich) ;
285  n_evol.update(nn_file, j, the_time[j]) ;
286  }
287  }
288 
289  // beta
290  for (int j=jmin; j<=jtime; j++) {
291  fread_be(&indicator, sizeof(int), 1, fich) ;
292  if (indicator == 1) {
293  Vector beta_file(mp, triad, fich) ;
294  beta_evol.update(beta_file, j, the_time[j]) ;
295  }
296  }
297 
298  // Case of a full reading
299  // -----------------------
300  if (!partial_read) {
301  cout <<
302  "Time_slice constructor from file: the case of full reading\n"
303  << " is not ready yet !" << endl ;
304  abort() ;
305  }
306 
307  set_der_0x0() ;
308 
309 }
310 
311 
312 // Copy constructor
313 // ----------------
315  : depth(tin.depth),
316  scheme_order(tin.scheme_order),
317  jtime(tin.jtime),
318  the_time(tin.the_time),
319  gam_dd_evol(tin.gam_dd_evol),
320  gam_uu_evol(tin.gam_uu_evol),
321  k_dd_evol(tin.k_dd_evol),
322  k_uu_evol(tin.k_uu_evol),
323  n_evol(tin.n_evol),
324  beta_evol(tin.beta_evol),
325  trk_evol(tin.trk_evol),
326  adm_mass_evol(tin.adm_mass_evol) {
327 
328  set_der_0x0() ;
329 }
330 
331 // Special constructor for derived classes
332 //----------------------------------------
334  : depth(depth_in),
335  scheme_order(depth_in-1),
336  jtime(0),
337  the_time(0., depth_in),
338  gam_dd_evol(depth_in),
339  gam_uu_evol(depth_in),
340  k_dd_evol(depth_in),
341  k_uu_evol(depth_in),
342  n_evol(depth_in),
343  beta_evol(depth_in),
344  trk_evol(depth_in),
345  adm_mass_evol() {
346 
347  set_der_0x0() ;
348 }
349 
350 
351 
352 
353  //--------------//
354  // Destructor //
355  //--------------//
356 
358 
360 
361 }
362 
363  //---------------------//
364  // Memory management //
365  //---------------------//
366 
367 void Time_slice::del_deriv() const {
368 
369  if (p_gamma != 0x0) delete p_gamma ;
370 
371  set_der_0x0() ;
372 
373  adm_mass_evol.downdate(jtime) ;
374 }
375 
376 
378 
379  p_gamma = 0x0 ;
380 
381 }
382 
383 
384  //-----------------------//
385  // Mutators / assignment //
386  //-----------------------//
387 
389 
390  depth = tin.depth;
391  scheme_order = tin.scheme_order ;
392  jtime = tin.jtime;
393  the_time = tin.the_time;
394  gam_dd_evol = tin.gam_dd_evol;
395  gam_uu_evol = tin.gam_uu_evol;
396  k_dd_evol = tin.k_dd_evol;
397  k_uu_evol = tin.k_uu_evol;
398  n_evol = tin.n_evol;
399  beta_evol = tin.beta_evol;
400  trk_evol = tin.trk_evol ;
401 
402  del_deriv() ;
403 
404 }
405 
406 
407  //------------------//
408  // output //
409  //------------------//
410 
411 ostream& Time_slice::operator>>(ostream& flux) const {
412 
413  flux << "\n------------------------------------------------------------\n"
414  << "Lorene class : " << typeid(*this).name() << '\n' ;
415  flux << "Number of stored slices : " << depth
416  << " order of time scheme : " << scheme_order << '\n'
417  << "Time label t = " << the_time[jtime]
418  << " index of time step j = " << jtime << '\n' << '\n' ;
419  if (adm_mass_evol.is_known(jtime)) {
420  flux << "ADM mass : " << adm_mass() << endl ;
421  }
422 
423  flux << "Max. of absolute values of the various fields in each domain: \n" ;
424  if (gam_dd_evol.is_known(jtime)) {
425  maxabs( gam_dd_evol[jtime], "gam_{ij}", flux) ;
426  }
427  if (gam_uu_evol.is_known(jtime)) {
428  maxabs( gam_uu_evol[jtime], "gam^{ij}", flux) ;
429  }
430  if (k_dd_evol.is_known(jtime)) {
431  maxabs( k_dd_evol[jtime], "K_{ij}", flux) ;
432  }
433  if (k_uu_evol.is_known(jtime)) {
434  maxabs( k_uu_evol[jtime], "K^{ij}", flux) ;
435  }
436  if (n_evol.is_known(jtime)) {
437  maxabs( n_evol[jtime], "N", flux) ;
438  }
439  if (beta_evol.is_known(jtime)) {
440  maxabs( beta_evol[jtime], "beta^i", flux) ;
441  }
442  if (trk_evol.is_known(jtime)) {
443  maxabs( trk_evol[jtime], "tr K", flux) ;
444  }
445 
446  if (p_gamma != 0x0) flux << "Metric gamma is up to date" << endl ;
447 
448  return flux ;
449 
450 }
451 
452 
453 ostream& operator<<(ostream& flux, const Time_slice& sigma) {
454 
455  sigma >> flux ;
456  return flux ;
457 
458 }
459 
460 
461 void Time_slice::save(const char* rootname) const {
462 
463  // Opening of file
464  // ---------------
465  char* filename = new char[ strlen(rootname)+10 ] ;
466  strcpy(filename, rootname) ;
467  char nomj[7] ;
468  sprintf(nomj, "%06d", jtime) ;
469  strcat(filename, nomj) ;
470  strcat(filename, ".d") ;
471 
472  FILE* fich = fopen(filename, "w") ;
473  if (fich == 0x0) {
474  cout << "Problem in opening file " << filename << " ! " << endl ;
475  perror(" reason") ;
476  abort() ;
477  }
478 
479  // Write grid, mapping, triad and depth
480  // ------------------------------------
481  const Map& map = nn().get_mp() ;
482  const Mg3d& mgrid = *(map.get_mg()) ;
483  const Base_vect& triad = *(beta().get_triad()) ;
484 
485  mgrid.sauve(fich) ;
486  map.sauve(fich) ;
487  triad.sauve(fich) ;
488 
489  fwrite_be(&depth, sizeof(int), 1, fich) ;
490 
491  // Write all binary data by means of virtual function sauve
492  // --------------------------------------------------------
493  bool partial_save = false ;
494  sauve(fich, partial_save) ;
495 
496  // Close the file
497  // --------------
498 
499  fclose(fich) ;
500 
501  delete [] filename ;
502 
503 }
504 
505 
506 
507 void Time_slice::sauve(FILE* fich, bool partial_save) const {
508 
509  // Writing various integer parameters
510  // ----------------------------------
511 
512  fwrite_be(&depth, sizeof(int), 1, fich) ;
513  fwrite_be(&scheme_order, sizeof(int), 1, fich) ;
514  fwrite_be(&jtime, sizeof(int), 1, fich) ;
515 
516  // Writing the_time
517  // ----------------
518  int jmin = jtime - depth + 1 ;
519  for (int j=jmin; j<=jtime; j++) {
520  int indicator = (the_time.is_known(j)) ? 1 : 0 ;
521  fwrite_be(&indicator, sizeof(int), 1, fich) ;
522  if (indicator == 1) {
523  double xx = the_time[j] ;
524  fwrite_be(&xx, sizeof(double), 1, fich) ;
525  }
526  }
527 
528  // Writing of various fields
529  // -------------------------
530 
531  // N
532  nn() ; // forces the update at the current time step
533  for (int j=jmin; j<=jtime; j++) {
534  int indicator = (n_evol.is_known(j)) ? 1 : 0 ;
535  fwrite_be(&indicator, sizeof(int), 1, fich) ;
536  if (indicator == 1) n_evol[j].sauve(fich) ;
537  }
538 
539  // beta
540  beta() ; // forces the update at the current time step
541  for (int j=jmin; j<=jtime; j++) {
542  int indicator = (beta_evol.is_known(j)) ? 1 : 0 ;
543  fwrite_be(&indicator, sizeof(int), 1, fich) ;
544  if (indicator == 1) beta_evol[j].sauve(fich) ;
545  }
546 
547  // Case of a complete save
548  // -----------------------
549  if (!partial_save) {
550 
551  cout << "Time_slice::sauve: the full writing is not ready yet !"
552  << endl ;
553  abort() ;
554  }
555 
556 
557 }
558 
559 
560 
561 
562 
563 
564 
565 
566 
567 
568 
569 
570 
571 
572 }
Cartesian vectorial bases (triads).
Definition: base_vect.h:201
Spherical orthonormal vectorial bases (triads).
Definition: base_vect.h:308
Vectorial bases (triads) with respect to which the tensorial components are defined.
Definition: base_vect.h:105
virtual void sauve(FILE *) const
Save in a file.
Definition: base_vect.C:150
Time evolution with partial storage (*** under development ***).
Definition: evolution.h:371
virtual void update(const TyT &new_value, int j, double time_j)
Sets a new value at a given time step.
bool is_known(int j) const
Checks whether the value a given time step has been set.
Definition: evolution.C:335
Base class for coordinate mappings.
Definition: map.h:670
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:765
const Metric_flat & flat_met_cart() const
Returns the flat metric associated with the Cartesian coordinates and with components expressed in th...
Definition: map.C:331
virtual void sauve(FILE *) const
Save in a file.
Definition: map.C:224
const Metric_flat & flat_met_spher() const
Returns the flat metric associated with the spherical coordinates and with components expressed in th...
Definition: map.C:321
virtual const Sym_tensor & cov() const
Read-only access to the covariant representation.
Definition: metric_flat.C:134
Metric for tensor calculation.
Definition: metric.h:90
Multi-domain grid.
Definition: grilles.h:273
void sauve(FILE *fd, bool save_base=false) const
Saves into a file.
Definition: mg3d.C:371
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:387
void set_etat_one()
Sets the logical state to ETATUN (one).
Definition: scalar.C:334
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field.
Definition: scalar.C:784
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: scalar.C:324
Class intended to describe valence-2 symmetric tensors.
Definition: sym_tensor.h:223
Spacelike time slice of a 3+1 spacetime.
Definition: time_slice.h:173
int jtime
Time step index of the latest slice.
Definition: time_slice.h:190
Evolution_std< Scalar > trk_evol
Values at successive time steps of the trace K of the extrinsic curvature.
Definition: time_slice.h:224
Metric * p_gamma
Pointer on the induced metric at the current time step (jtime)
Definition: time_slice.h:239
Evolution_std< Sym_tensor > k_dd_evol
Values at successive time steps of the covariant components of the extrinsic curvature tensor .
Definition: time_slice.h:208
virtual void sauve(FILE *fich, bool partial_save) const
Total or partial saves in a binary file.
Definition: time_slice.C:507
void operator=(const Time_slice &)
Assignment to another Time_slice.
Definition: time_slice.C:388
Evolution_full< Tbl > adm_mass_evol
ADM mass at each time step, since the creation of the slice.
Definition: time_slice.h:233
virtual const Scalar & nn() const
Lapse function N at the current time step (jtime )
void save(const char *rootname) const
Saves in a binary file.
Definition: time_slice.C:461
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: time_slice.C:377
Evolution_std< Sym_tensor > gam_dd_evol
Values at successive time steps of the covariant components of the induced metric .
Definition: time_slice.h:198
Evolution_std< Sym_tensor > gam_uu_evol
Values at successive time steps of the contravariant components of the induced metric .
Definition: time_slice.h:203
int scheme_order
Order of the finite-differences scheme for the computation of time derivatives.
Definition: time_slice.h:187
virtual const Vector & beta() const
shift vector at the current time step (jtime )
virtual ~Time_slice()
Destructor.
Definition: time_slice.C:357
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: time_slice.C:367
int depth
Number of stored time slices.
Definition: time_slice.h:179
Evolution_std< Scalar > n_evol
Values at successive time steps of the lapse function N.
Definition: time_slice.h:216
virtual double adm_mass() const
Returns the ADM mass (geometrical units) at the current step.
Time_slice(const Scalar &lapse_in, const Vector &shift_in, const Sym_tensor &gamma_in, const Sym_tensor &kk_in, int depth_in=3)
General constructor (Hamiltonian-like).
Definition: time_slice.C:112
Evolution_std< double > the_time
Time label of each slice.
Definition: time_slice.h:193
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator<<).
Definition: time_slice.C:411
Evolution_std< Sym_tensor > k_uu_evol
Values at successive time steps of the contravariant components of the extrinsic curvature tensor .
Definition: time_slice.h:213
Evolution_std< Vector > beta_evol
Values at successive time steps of the shift vector .
Definition: time_slice.h:219
Tensor field of valence 1.
Definition: vector.h:188
int fread_be(int *aa, int size, int nb, FILE *fich)
Reads integer(s) from a binary file according to the big endian convention.
Definition: fread_be.C:69
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.
Definition: fwrite_be.C:70
const Map & get_mp() const
Returns the mapping.
Definition: tensor.h:861
int get_index_type(int i) const
Gives the type (covariant or contravariant) of the index number i .
Definition: tensor.h:886
Tensor trace(int ind1, int ind2) const
Trace on two different type indices.
virtual void set_etat_zero()
Sets the logical state of all components to ETATZERO (zero state).
Definition: tensor.C:497
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Definition: tensor.h:866
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.
Lorene prototypes.
Definition: app_hor.h:64