MBDyn-1.7.3
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups
varstep.cc
Go to the documentation of this file.
1 /* $Header: /var/cvs/mbdyn/mbdyn/mbdyn-1.0/mbdyn/base/varstep.cc,v 1.15 2017/01/12 14:46:11 masarati Exp $ */
2 /*
3  * MBDyn (C) is a multibody analysis code.
4  * http://www.mbdyn.org
5  *
6  * Copyright (C) 1996-2017
7  *
8  * Pierangelo Masarati <masarati@aero.polimi.it>
9  * Paolo Mantegazza <mantegazza@aero.polimi.it>
10  *
11  * Dipartimento di Ingegneria Aerospaziale - Politecnico di Milano
12  * via La Masa, 34 - 20156 Milano, Italy
13  * http://www.aero.polimi.it
14  *
15  * Changing this copyright notice is forbidden.
16  *
17  * This program is free software; you can redistribute it and/or modify
18  * it under the terms of the GNU General Public License as published by
19  * the Free Software Foundation (version 2 of the License).
20  *
21  *
22  * This program is distributed in the hope that it will be useful,
23  * but WITHOUT ANY WARRANTY; without even the implied warranty of
24  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
25  * GNU General Public License for more details.
26  *
27  * You should have received a copy of the GNU General Public License
28  * along with this program; if not, write to the Free Software
29  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
30  */
31 
32 /* fixed step file driver */
33 
34 #include "mbconfig.h" /* This goes first in every *.c,*.cc file */
35 
36 #include <fstream>
37 
38 #include "dataman.h"
39 #include "filedrv.h"
40 #include "varstep.h"
41 #include "solver.h"
42 #include "bisec.h"
43 
44 /* VariableStepFileDrive - begin */
45 
46 static const std::vector<doublereal> v0;
47 
49  const DriveHandler* pDH,
50  const char* const sFileName,
51  integer ind, bool bl, bool pz, Drive::Bailout bo)
52 : FileDrive(uL, pDH, sFileName, ind, v0),
53 iNumSteps(-1), iCurrStep(-1),
54 bLinear(bl), bPadZeroes(pz), boWhen(bo), pd(0), pvd(0)
55 {
56  ASSERT(iNumDrives > 0);
57  ASSERT(sFileName != NULL);
58 
59  std::ifstream in(sFileName);
60  if (!in) {
61  silent_cerr("can't open file \""
62  << sFileName << "\"" << std::endl);
64  }
65 
66  /*
67  * Mangia gli eventuali commenti iniziali
68  */
69  char c = '\0';
70  while (in.get(c), c == '#') {
71  char buf[1024];
72 
73  do {
74  in.getline(buf, sizeof(buf));
75  } while (strlen(buf) == STRLENOF(buf) && buf[STRLENOF(buf)] != '\n');
76  }
77 
78  if (c != '#') {
79  in.putback(c);
80  }
81 
82  if (iNumSteps == -1) {
83  std::streampos pos = in.tellg();
84 
85  integer ins;
86  for (ins = 0; !in.eof(); ins++) {
87  char buf[1024];
88 
89  do {
90  in.getline(buf, sizeof(buf));
91  } while (strlen(buf) == STRLENOF(buf) && buf[STRLENOF(buf)] != '\n');
92  }
93  iNumSteps = --ins;
94 
95  in.clear();
96  in.seekg(pos);
97 
98  silent_cout("VariableStepFileDrive(" << uL << "): "
99  "counted " << ins << " steps" << std::endl);
100  }
101 
104 
105  for (integer i = iNumDrives + 1; i-- > 0; ) {
106  pvd[i] = pd + i*iNumSteps;
107  }
108 
109  for (integer j = 0; j < iNumSteps; j++) {
110  // 0 -> iNumDrives to account for time
111  for (integer i = 0; i <= iNumDrives; i++) {
112  in >> pvd[i][j];
113  if (in.eof()) {
114  silent_cerr("unexpected end of file '"
115  << sFileName << '\'' << std::endl);
117  }
118 
119  int c;
120  for (c = in.get(); isspace(c); c = in.get()) {
121  if (c == '\n') {
122  if (i == 0) {
123  silent_cerr("unexpected end of line #" << j + 1 << " after time=" << pvd[0][j] << ", column #" << i + 1 << " of file '"
124  << sFileName << '\'' << std::endl);
126  }
127  if (i != iNumDrives) {
128  silent_cerr("unexpected end of line #" << j + 1 << ", time=" << pvd[0][j] << " after channel #" << i << ", column #" << i + 1 << " of file '"
129  << sFileName << '\'' << std::endl);
131  }
132  break;
133  }
134  }
135 
136  if (i == iNumDrives && c != '\n') {
137  silent_cerr("missing end-of-line at line #" << j + 1 << ", time=" << pvd[0][j] << " of file '"
138  << sFileName << '\'' << std::endl);
140  }
141 
142  in.putback(c);
143  }
144 
145  if (j > 0) {
146  if (pvd[0][j] <= pvd[0][j - 1]) {
147  silent_cerr("time[" << j << "]=" << pvd[0][j]
148  << " <= time[" << j - 1 << "]=" << pvd[0][j - 1]
149  << " in file '" << sFileName << "'" << std::endl);
151  }
152  }
153  }
154 
155  // All data is available, so initialize the buffer accordingly
156  // use bisection to initialize iCurrStep
157  doublereal dTime = pDH->dGetTime();
158  iCurrStep = bisec(pvd[0], dTime, 0, iNumSteps - 1);
159  if (iCurrStep < 0) {
160  iCurrStep++;
161  }
162 
163  ServePending(dTime);
164 }
165 
167 {
168  SAFEDELETEARR(pd);
170 }
171 
172 
173 /* Scrive il contributo del DriveCaller al file di restart */
174 std::ostream&
175 VariableStepFileDrive::Restart(std::ostream& out) const
176 {
177  return out << "0. /* VariableStepFileDrive: not implemented yet! */"
178  << std::endl;
179 }
180 
181 void
183 {
184  if (t <= pvd[0][0]) {
185  if (t < pvd[0][0] && (boWhen & Drive::BO_LOWER)) {
186  throw Solver::EndOfSimulation(EXIT_SUCCESS,
188  "A variable step file drive lower bound is halting the simulation");
189  }
190 
191  if (bPadZeroes) {
192  for (int i = 1; i <= iNumDrives; i++) {
193  pdVal[i] = 0.;
194  }
195 
196  } else {
197  for (int i = 1; i <= iNumDrives; i++) {
198  pdVal[i] = pvd[i][0];
199  }
200  }
201 
202  } else if (t >= pvd[0][iNumSteps - 1]) {
203  if (t > pvd[0][iNumSteps - 1] && (boWhen & Drive::BO_UPPER)) {
204  throw Solver::EndOfSimulation(EXIT_SUCCESS,
206  "A variable step file drive upper bound is halting the simulation");
207  }
208 
209  if (bPadZeroes) {
210  for (int i = 1; i <= iNumDrives; i++) {
211  pdVal[i] = 0.;
212  }
213 
214  } else {
215  for (int i = 1; i <= iNumDrives; i++) {
216  pdVal[i] = pvd[i][iNumSteps - 1];
217  }
218  }
219 
220  } else {
221  // look for step exactly before
222  // note: use linear search, under the assumption
223  // that we always start from a relatively close guess
224  while (pvd[0][iCurrStep] > t) {
225  iCurrStep--;
226  }
227 
228  // no need to check for under/uverflow
229  ASSERT(iCurrStep >= 0);
230 
231  while (pvd[0][iCurrStep + 1] <= t) {
232  iCurrStep++;
233  }
235 
236  integer j1 = iCurrStep;
237  if (bLinear) {
238  integer j2 = j1 + 1;
239  doublereal dt1 = pvd[0][j1];
240  doublereal dt2 = pvd[0][j2];
241  doublereal dw1 = (dt2 - t)/(dt2 - dt1);
242  doublereal dw2 = (t - dt1)/(dt2 - dt1);
243 
244  for (int i = 1; i <= iNumDrives; i++) {
245  pdVal[i] = pvd[i][j2]*dw2 + pvd[i][j1]*dw1;
246  }
247 
248  } else {
249  for (int i = 1; i <= iNumDrives; i++) {
250  pdVal[i] = pvd[i][j1];
251  }
252  }
253  }
254 }
255 
256 /* VariableStepFileDrive - end */
257 
258 
259 /* legge i drivers tipo fixed step file */
260 
261 Drive *
262 VariableStepDR::Read(unsigned uLabel, const DataManager *pDM, MBDynParser& HP)
263 {
264  integer idrives = HP.GetInt();
265  if (idrives <= 0) {
266  silent_cerr("VariableStepFileDrive(" << uLabel << "): "
267  "invalid channels number " << idrives
268  << " at line " << HP.GetLineData()
269  << std::endl);
271  }
272 
273  bool bl(true);
274  if (HP.IsKeyWord("interpolation")) {
275  if (HP.IsKeyWord("const")) {
276  bl = false;
277 
278  } else if (!HP.IsKeyWord("linear")) {
279  silent_cerr("VariableStepFileDrive(" << uLabel << "): "
280  "unknown value for \"interpolation\" "
281  "at line " << HP.GetLineData() << std::endl);
283  }
284  }
285 
286  bool pz(true);
288 
289  if (HP.IsKeyWord("pad" "zeros") || HP.IsKeyWord("pad" "zeroes")) {
290  if (!HP.GetYesNo(pz)) {
291  silent_cerr("VariableStepFileDrive(" << uLabel << "): "
292  "unknown value for \"pad zeros\" "
293  "at line " << HP.GetLineData() << std::endl);
295  }
296 
297  } else if (HP.IsKeyWord("bailout")) {
298  if (HP.IsKeyWord("none")) {
299  bo = Drive::BO_NONE;
300 
301  } else if (HP.IsKeyWord("upper")) {
302  bo = Drive::BO_UPPER;
303 
304  } else if (HP.IsKeyWord("lower")) {
305  bo = Drive::BO_LOWER;
306 
307  } else if (HP.IsKeyWord("any")) {
308  bo = Drive::BO_ANY;
309 
310  } else {
311  silent_cerr("VariableStepFileDrive(" << uLabel << "): "
312  "invalid bailout parameter "
313  "at line " << HP.GetLineData()
314  << std::endl);
316  }
317  }
318 
319  const char* filename = HP.GetFileName();
320 
321  Drive* pDr = NULL;
324  VariableStepFileDrive(uLabel, pDM->pGetDrvHdl(),
325  filename, idrives, bl, pz, bo));
326 
327  return pDr;
328 }
329 
static char * filename
Definition: ann_sf.c:71
#define MBDYN_EXCEPT_ARGS
Definition: except.h:63
virtual integer GetInt(integer iDefval=0)
Definition: parser.cc:1050
Definition: drive.h:89
#define SAFEDELETEARR(pnt)
Definition: mynewmem.h:713
virtual const char * GetFileName(enum Delims Del=DEFAULTDELIM)
Definition: parsinc.cc:673
Bailout
Definition: drive.h:101
const DriveHandler * pGetDrvHdl(void) const
Definition: dataman.h:340
virtual Drive * Read(unsigned uLabel, const DataManager *pDM, MBDynParser &HP)
Definition: varstep.cc:262
integer iNumDrives
Definition: filedrv.h:47
virtual void ServePending(const doublereal &t)
Definition: varstep.cc:182
doublereal * pd
Definition: varstep.h:49
doublereal * pdVal
Definition: filedrv.h:48
static const std::vector< doublereal > v0
Definition: varstep.cc:46
virtual bool IsKeyWord(const char *sKeyWord)
Definition: parser.cc:910
VariableStepFileDrive(unsigned int uL, const DriveHandler *pDH, const char *const sFileName, integer id, bool bl, bool pz, Drive::Bailout bo)
Definition: varstep.cc:48
#define ASSERT(expression)
Definition: colamd.c:977
#define SAFENEWWITHCONSTRUCTOR(pnt, item, constructor)
Definition: mynewmem.h:698
static std::stack< cleanup * > c
Definition: cleanup.cc:59
#define STRLENOF(s)
Definition: mbdyn.h:166
virtual bool GetYesNo(bool &bRet)
Definition: parser.cc:1022
doublereal ** pvd
Definition: varstep.h:50
doublereal dGetTime(void) const
Definition: drive.h:386
#define SAFENEWARR(pnt, item, sz)
Definition: mynewmem.h:701
virtual ~VariableStepFileDrive(void)
Definition: varstep.cc:166
static doublereal buf[BUFSIZE]
Definition: discctrl.cc:333
double doublereal
Definition: colamd.c:52
long int integer
Definition: colamd.c:51
virtual HighParser::ErrOut GetLineData(void) const
Definition: parsinc.cc:697
virtual std::ostream & Restart(std::ostream &out) const
Definition: varstep.cc:175