LORENE
misc.C
1 /*
2  * Various utilities...
3  *
4  */
5 
6 /*
7  * Copyright (c) 2018 Jerome Novak
8  *
9  * This file is part of LORENE.
10  *
11  * LORENE is free software; you can redistribute it and/or modify
12  * it under the terms of the GNU General Public License version 2
13  * as published by the Free Software Foundation.
14  *
15  * LORENE is distributed in the hope that it will be useful,
16  * but WITHOUT ANY WARRANTY; without even the implied warranty of
17  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18  * GNU General Public License for more details.
19  *
20  * You should have received a copy of the GNU General Public License
21  * along with LORENE; if not, write to the Free Software
22  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
23  *
24  */
25 
26 char misc_C[] = "$Header: /cvsroot/Lorene/C++/Source/Non_class_members/Utilities/misc.C,v 1.2 2019/12/02 14:51:37 j_novak Exp $" ;
27 
28 /*
29  * $Id: misc.C,v 1.2 2019/12/02 14:51:37 j_novak Exp $
30  * $Log: misc.C,v $
31  * Revision 1.2 2019/12/02 14:51:37 j_novak
32  * Moved piecewise parabolic computation of dpdnb to a separate function.
33  *
34  * Revision 1.1 2018/12/05 15:03:20 j_novak
35  * New Mg3d constructor from a formatted file.
36  *
37  * Revision 1.4 2003/10/19 20:01:10 e_gourgoulhon
38  * Template file
39  *
40  * $Header: /cvsroot/Lorene/C++/Source/Non_class_members/Utilities/misc.C,v 1.2 2019/12/02 14:51:37 j_novak Exp $
41  *
42  */
43 
44 // Lorene headers
45 #include "tbl.h"
46 
47 
48 namespace Lorene {
49 
50  // Searches the file 'infile' for a given 'pattern'. The file stream is
51  // positionned just after the occurrence of the pattern.
52  //=====================================================================
53  bool search_file(ifstream& infile, const string& pattern) {
54  string line ;
55  while( getline(infile, line) ) {
56  if (line.find(pattern, 0) != string::npos)
57  return true ;
58  }
59  return false ;
60  }
61 
62  // Computes the first derivative of the 2nd Tbl, with respect to the 1st one
63  // using piecewise parabolic scheme.
64  void compute_derivative(const Tbl& xx, const Tbl& yy, Tbl& dydx)
65  {
66  double y0, y1, y2, x0, x1, x2, der;
67  int np = xx.get_dim(0) ;
68  assert ((yy.get_dim(0) == np) && (dydx.get_dim(0) == np)) ;
69 
70  // special case: i=0
71 
72  y0 = yy(0);
73  y1 = yy(1);
74  y2 = yy(2);
75 
76  x0 = xx(0);
77  x1 = xx(1);
78  x2 = xx(2);
79 
80  der = y0*(2*x0-x1-x2)/(x0-x1)/(x0-x2) +
81  y1*(x0-x2)/(x1-x0)/(x1-x2) +
82  y2*(x0-x1)/(x2-x0)/(x2-x1) ;
83 
84  dydx.set(0) = der ;
85 
86  for(int i=1;i<np-1;i++) {
87 
88  y0 = yy(i-1);
89  y1 = yy(i);
90  y2 = yy(i+1);
91 
92  x0 = xx(i-1);
93  x1 = xx(i);
94  x2 = xx(i+1);
95 
96  der = y0*(x1-x2)/(x0-x1)/(x0-x2) +
97  y1*(2*x1-x0-x2)/(x1-x0)/(x1-x2) +
98  y2*(x1-x0)/(x2-x0)/(x2-x1) ;
99 
100  dydx.set(i) = der ;
101 
102  }
103 
104  // special case: i=np-1
105 
106  y0 = yy(np-3);
107  y1 = yy(np-2);
108  y2 = yy(np-1);
109 
110  x0 = xx(np-3);
111  x1 = xx(np-2);
112  x2 = xx(np-1);
113 
114  der = y0*(x2-x1)/(x0-x1)/(x0-x2) +
115  y1*(x2-x0)/(x1-x0)/(x1-x2) +
116  y2*(2*x2-x0-x1)/(x2-x0)/(x2-x1) ;
117 
118  dydx.set(np-1) = der ;
119 
120  return ;
121  }
122 
123 } // End of namespace Lorene
Lorene prototypes.
Definition: app_hor.h:67
double & set(int i)
Read/write of a particular element (index i) (1D case)
Definition: tbl.h:301
bool search_file(ifstream &infile, const string &pattern)
A function that searches for a pattern in a file and places the file stream after the found pattern...
Definition: misc.C:53
void compute_derivative(const Tbl &xx, const Tbl &ff, Tbl &dfdx)
Derives a function defined on an unequally-spaced grid, approximating it by piecewise parabolae...
Definition: misc.C:64
int get_dim(int i) const
Gives the i-th dimension (ie dim.dim[i])
Definition: tbl.h:423
Basic array class.
Definition: tbl.h:164