MBDyn-1.7.3
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups
c81merge.cc
Go to the documentation of this file.
1 /* $Header: /var/cvs/mbdyn/mbdyn/mbdyn-1.0/utils/c81merge.cc,v 1.27 2017/01/12 15:10:27 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  *
10  * Dipartimento di Ingegneria Aerospaziale - Politecnico di Milano
11  * via La Masa, 34 - 20156 Milano, Italy
12  * http://www.aero.polimi.it
13  *
14  * Changing this copyright notice is forbidden.
15  *
16  * This program is free software; you can redistribute it and/or modify
17  * it under the terms of the GNU General Public License as published by
18  * the Free Software Foundation (version 2 of the License).
19  *
20  *
21  * This program is distributed in the hope that it will be useful,
22  * but WITHOUT ANY WARRANTY; without even the implied warranty of
23  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
24  * GNU General Public License for more details.
25  *
26  * You should have received a copy of the GNU General Public License
27  * along with this program; if not, write to the Free Software
28  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
29  */
30 
31 #include "mbconfig.h"
32 
33 #include <cstdlib>
34 #include <unistd.h>
35 
36 #include <iostream>
37 #include <fstream>
38 #include "ac/getopt.h"
39 
40 #include "myassert.h"
41 #include "mynewmem.h"
42 
43 #include "aerodata.h"
44 #include "aerodc81.h"
45 #include "c81data.h"
46 
47 static void
49 {
50  for (int i = 0; i < n; i++) {
51  v[i] = w1*v1[i] + w2*v2[i];
52  }
53 }
54 
55 static void
56 merge_data_from(int NA, int NM, int NM_to,
57  doublereal wfrom, doublereal wto, doublereal *workv,
58  doublereal *m_from, doublereal *m_to, doublereal *m_dst,
59  doublereal *a_from, doublereal *a_to, doublereal *a_dst)
60 {
61  int mt = 0;
62 
63  for (int m = 0; m < NM; m++) {
64  m_dst[m] = m_from[m];
65  }
66 
67  for (int a = 0; a < NA; a++) {
68  a_dst[a] = a_from[a];
69  }
70 
71  for (int m = 0; m < NM; m++) {
72  for (; mt < NM_to; mt++) {
73  if (m_to[mt] >= m_dst[m]) {
74  break;
75  }
76  }
77 
78  /* last, step back and maintain */
79  if (mt == NM_to) {
80  mt = NM_to - 1;
81  merge_vecs(NA, &a_from[NA*(1 + m)], wfrom,
82  &a_to[NA*(1 + mt)], wto,
83  &a_dst[NA*(1 + m)]);
84 
85  /* first, or exact match */
86  } else if (mt == 0 || m_to[mt] == m_dst[m]) {
87  merge_vecs(NA, &a_from[NA*(1 + m)], wfrom,
88  &a_to[NA*(1 + mt)], wto,
89  &a_dst[NA*(1 + m)]);
90 
91  /* need interpolation */
92  } else {
93  doublereal dd = m_to[mt] - m_to[mt - 1];
94  doublereal d1 = (m_to[mt] - m_dst[m])/dd;
95  doublereal d2 = (m_dst[m] - m_to[mt - 1])/dd;
96  merge_vecs(NA, &a_to[NA*(mt)], d1,
97  &a_to[NA*(mt + 1)], d2,
98  workv);
99  merge_vecs(NA, &a_from[NA*(1 + m)], wfrom,
100  workv, wto,
101  &a_dst[NA*(1 + m)]);
102  }
103  }
104 }
105 
106 static int
107 c81_data_merge_from(c81_data *data_from, doublereal wfrom, c81_data *data_to, doublereal wto, c81_data *data)
108 {
109  *data = *data_from;
110 
111  data->ml = new doublereal[data->NML];
112  data->al = new doublereal[data->NAL*(data->NML + 1)];
113 
114  data->md = new doublereal[data->NMD];
115  data->ad = new doublereal[data->NAD*(data->NMD + 1)];
116 
117  data->mm = new doublereal[data->NMM];
118  data->am = new doublereal[data->NAM*(data->NMM + 2)];
119 
120  int NA, NM, NM_to;
121  doublereal *workv = &data->am[data->NAM*(data->NMM + 1)],
122  *a_from, *a_to, *a_dst, *m_from, *m_to, *m_dst;
123 
124  a_from = data_from->al;
125  a_to = data_to->al;
126  a_dst = data->al;
127  m_from = data_from->ml;
128  m_to = data_to->ml;
129  m_dst = data->ml;
130  NA = data->NAL;
131  NM = data->NML;
132  NM_to = data_to->NML;
133  merge_data_from(NA, NM, NM_to, wfrom, wto, workv, m_from, m_to, m_dst, a_from, a_to, a_dst);
134 
135  a_from = data_from->ad;
136  a_to = data_to->ad;
137  a_dst = data->ad;
138  m_from = data_from->md;
139  m_to = data_to->md;
140  m_dst = data->md;
141  NA = data->NAD;
142  NM = data->NMD;
143  NM_to = data_to->NMD;
144  merge_data_from(NA, NM, NM_to, wfrom, wto, workv, m_from, m_to, m_dst, a_from, a_to, a_dst);
145 
146  a_from = data_from->am;
147  a_to = data_to->am;
148  a_dst = data->am;
149  m_from = data_from->mm;
150  m_to = data_to->mm;
151  m_dst = data->mm;
152  NA = data->NAM;
153  NM = data->NMM;
154  NM_to = data_to->NMM;
155  merge_data_from(NA, NM, NM_to, wfrom, wto, workv, m_from, m_to, m_dst, a_from, a_to, a_dst);
156 
157  // c81_do_data_stall(data);
158 
159  return 0;
160 }
161 
162 static void
163 merge_data_to(int NA, int NM, int NM_from,
164  doublereal wfrom, doublereal wto, doublereal *workv,
165  doublereal *m_from, doublereal *m_to, doublereal *m_dst,
166  doublereal *a_from, doublereal *a_to, doublereal *a_dst)
167 {
168  int mf = 0;
169 
170  for (int m = 0; m < NM; m++) {
171  m_dst[m] = m_to[m];
172  }
173 
174  for (int a = 0; a < NA; a++) {
175  a_dst[a] = a_to[a];
176  }
177 
178  for (int m = 0; m < NM; m++) {
179  for (; mf < NM_from; mf++) {
180  if (m_from[mf] >= m_dst[m]) {
181  break;
182  }
183  }
184 
185  /* last, step back and maintain */
186  if (mf == NM_from) {
187  mf = NM_from - 1;
188  merge_vecs(NA, &a_from[NA*(1 + mf)], wfrom,
189  &a_to[NA*(1 + m)], wto,
190  &a_dst[NA*(1 + m)]);
191 
192  /* first, or exact match */
193  } else if (mf == 0 || m_from[mf] == m_dst[m]) {
194  merge_vecs(NA, &a_from[NA*(1 + mf)], wfrom,
195  &a_to[NA*(1 + m)], wto,
196  &a_dst[NA*(1 + m)]);
197 
198  /* need interpolation */
199  } else {
200  doublereal dd = m_from[mf] - m_from[mf - 1];
201  doublereal d1 = (m_from[mf] - m_dst[m])/dd;
202  doublereal d2 = (m_dst[m] - m_from[mf - 1])/dd;
203  merge_vecs(NA, &a_from[NA*(mf)], d1,
204  &a_from[NA*(mf + 1)], d2,
205  workv);
206  merge_vecs(NA, workv, wfrom,
207  &a_to[NA*(1 + m)], wto,
208  &a_dst[NA*(1 + m)]);
209  }
210  }
211 }
212 
213 static int
214 c81_data_merge_to(c81_data *data_from, doublereal wfrom, c81_data *data_to, doublereal wto, c81_data *data)
215 {
216  *data = *data_to;
217 
218  data->ml = new doublereal[data->NML];
219  data->al = new doublereal[data->NAL*(data->NML + 1)];
220 
221  data->md = new doublereal[data->NMD];
222  data->ad = new doublereal[data->NAD*(data->NMD + 1)];
223 
224  data->mm = new doublereal[data->NMM];
225  data->am = new doublereal[data->NAM*(data->NMM + 2)];
226 
227  int NA, NM, NM_from;
228  doublereal *workv = &data->am[data->NAM*(data->NMM + 1)],
229  *a_from, *a_to, *a_dst, *m_from, *m_to, *m_dst;
230 
231  a_from = data_from->al;
232  a_to = data_to->al;
233  a_dst = data->al;
234  m_from = data_from->ml;
235  m_to = data_to->ml;
236  m_dst = data->ml;
237  NA = data->NAL;
238  NM = data->NML;
239  NM_from = data_from->NML;
240  merge_data_to(NA, NM, NM_from, wfrom, wto, workv, m_from, m_to, m_dst, a_from, a_to, a_dst);
241 
242  a_from = data_from->ad;
243  a_to = data_to->ad;
244  a_dst = data->ad;
245  m_from = data_from->md;
246  m_to = data_to->md;
247  m_dst = data->md;
248  NA = data->NAD;
249  NM = data->NMD;
250  NM_from = data_from->NMD;
251  merge_data_to(NA, NM, NM_from, wfrom, wto, workv, m_from, m_to, m_dst, a_from, a_to, a_dst);
252 
253  a_from = data_from->am;
254  a_to = data_to->am;
255  a_dst = data->am;
256  m_from = data_from->mm;
257  m_to = data_to->mm;
258  m_dst = data->mm;
259  NA = data->NAM;
260  NM = data->NMM;
261  NM_from = data_from->NMM;
262  merge_data_to(NA, NM, NM_from, wfrom, wto, workv, m_from, m_to, m_dst, a_from, a_to, a_dst);
263 
264  // c81_do_data_stall(data);
265 
266  return 0;
267 }
268 
269 static int
270 c81_data_merge_both(c81_data *data_from, doublereal wfrom, c81_data *data_to, doublereal wto, c81_data *data)
271 {
272  return -1;
273 }
274 
275 /*
276  * uso:
277  */
278 int
279 main(int argc, char *argv[])
280 {
281  c81_data data_from,
282  data_to,
283  data;
284  char *name_from = 0,
285  *name_to = 0,
286  *name = 0,
287  header[31] = { '\0' },
288  *next;
289  doublereal dFrom = -1.,
290  dTo;
291  int rc = EXIT_SUCCESS;
292  doublereal tol = 1e-6;
293 
294  enum {
295  MODE_UNDEFINED, MODE_FROM, MODE_TO, MODE_BOTH
296  } mode = MODE_UNDEFINED;
297 
298  while (1) {
299  int opt = getopt(argc, argv, "f:hH:m:o:s:t:");
300 
301  if (opt == EOF) {
302  break;
303  }
304 
305  switch (opt) {
306  case 'f':
307  if (name_from != 0) {
308  silent_cerr("-f already set to \"" << name_from << "\"" << std::endl);
309  }
310  name_from = optarg;
311  break;
312 
313  case 'H':
314  if (header[0] != '\0') {
315  silent_cerr("-h already set to \"" << header << "\"" << std::endl);
316  }
317  memcpy(header, optarg, sizeof(header));
318  header[STRLENOF(header)] = '\0';
319  break;
320 
321  case 'm':
322  switch (mode) {
323  case MODE_FROM:
324  silent_cerr("-m already set to \"FROM\"" << std::endl);
325  break;
326 
327  case MODE_TO:
328  silent_cerr("-m already set to \"TO\"" << std::endl);
329  break;
330 
331  case MODE_BOTH:
332  silent_cerr("-m already set to \"BOTH\"" << std::endl);
333  break;
334 
335  case MODE_UNDEFINED:
336  break;
337  }
338 
339  if (strcasecmp(optarg, "from") == 0) {
340  mode = MODE_FROM;
341 
342  } else if (strcasecmp(optarg, "to") == 0) {
343  mode = MODE_TO;
344 
345  } else if (strcasecmp(optarg, "both") == 0) {
346  mode = MODE_BOTH;
347 
348  } else {
349  silent_cerr("unknown mode \"" << optarg << "\"" <<std::endl);
351  }
352  break;
353 
354  case 'o':
355  if (name != 0) {
356  silent_cerr("-o already set to \"" << name << "\"" << std::endl);
357  }
358  name = optarg;
359  break;
360 
361  case 's':
362  if (dFrom != -1.) {
363  silent_cerr("-s already set to " << dFrom << std::endl);
364  }
365  dFrom = strtod(optarg, &next);
366  if (next == optarg || next[0] != '\0') {
367  silent_cerr("illegal value \"" << optarg << "\" for -s option" << std::endl);
369  }
370  if (dFrom < 0. || dFrom > 1.) {
371  silent_cerr("-s " << optarg << " is out of bounds (0.,1.)" << std::endl);
373  }
374  break;
375 
376  case 't':
377  if (name_to != 0) {
378  silent_cerr("-t already set to \"" << name_to << "\"" << std::endl);
379  }
380  name_to = optarg;
381  break;
382 
383  default:
384  rc = EXIT_FAILURE;
385 
386  case 'h':
387  silent_cout(
388 "\n"
389 "c81merge: merges two c81 files in a given proportion\n"
390 "\n"
391 "\t-f first.c81\t" "first c81 data set\n"
392 "\t-h\t\t" "this message\n"
393 "\t-H header\t" "the header (max 30 chars)\n"
394 "\t-m {from|to|both}\t" "use alpha/mach of \"from\", \"to\" or both\n"
395 "\t-o out.c81\t" "output file name (stdout if missing)\n"
396 "\t-s s in [0,1]\t" "real number; result = s * first + (1 - s) * second\n"
397 "\t-t second.c81\t" "second c81 data set\n"
398 "\n"
399  );
400  exit(rc);
401  }
402  }
403 
404  if (name_from == 0) {
405  silent_cerr("missing required -f \"from\" parameter" << std::endl);
407  }
408 
409  if (name_to == 0) {
410  silent_cerr("missing required -t \"to\" parameter" << std::endl);
412  }
413 
414  if (dFrom == -1.) {
415  silent_cerr("missing required -s \"fraction\" parameter" << std::endl);
417  }
418 
419  if (mode == MODE_UNDEFINED) {
420  mode = MODE_BOTH;
421  }
422 
423  dTo = 1. - dFrom;
424 
425  std::ifstream in;
426 
427  /* from */
428  in.open(name_from);
429  if (!in) {
430  silent_cerr("unable to open file \"" << name_from << "\"" << std::endl);
432  }
433 
434  int ff_from = 0;
435  if (c81_data_read(in, &data_from, tol, &ff_from)) {
436  silent_cerr("unable to read c81 data from file \"" << name_from << "\"" << std::endl);
438  }
439  in.close();
440 
441  /* to */
442  in.open(name_to);
443  if (!in) {
444  silent_cerr("unable to open file \"" << name_to << "\"" << std::endl);
446  }
447 
448  int ff_to = 0;
449  if (c81_data_read(in, &data_to, tol, &ff_to)) {
450  silent_cerr("unable to read c81 data from file \"" << name_to << "\"" << std::endl);
452  }
453  in.close();
454 
455  /* consistency: at the moment, the two data sets must have exactly the same AoA pattern */
456  if (data_from.NAL != data_to.NAL) {
457  silent_cerr("number of AoA values for Cl differ ("
458  << data_from.NAL << " vs. " << data_to.NAL
459  << ")" << std::endl);
461  }
462 
463  for (int i = 0; i < data_from.NAL; i++) {
464  if (data_from.al[i] != data_to.al[i]) {
465  silent_cerr("AoA value " << i << " for Cl differs ("
466  << data_from.al[i] << " vs. "
467  << data_to.al[i] << ")"
468  << std::endl);
470  }
471  }
472 
473  if (data_from.NAD != data_to.NAD) {
474  silent_cerr("number of AoA values for Cd differ ("
475  << data_from.NAD << " vs. " << data_to.NAD
476  << ")" << std::endl);
478  }
479 
480  for (int i = 0; i < data_from.NAL; i++) {
481  if (data_from.ad[i] != data_to.ad[i]) {
482  silent_cerr("AoA value " << i << " for Cd differs ("
483  << data_from.ad[i] << " vs. "
484  << data_to.ad[i] << ")"
485  << std::endl);
487  }
488  }
489 
490  if (data_from.NAM != data_to.NAM) {
491  silent_cerr("number of AoA values for Cm differ ("
492  << data_from.NAM << " vs. " << data_to.NAM
493  << ")" << std::endl);
495  }
496 
497  for (int i = 0; i < data_from.NAL; i++) {
498  if (data_from.am[i] != data_to.am[i]) {
499  silent_cerr("AoA value " << i << " for Cm differs ("
500  << data_from.am[i] << " vs. "
501  << data_to.am[i] << ")"
502  << std::endl);
504  }
505  }
506 
507  /* merge */
508  switch (mode) {
509  case MODE_FROM:
510  rc = c81_data_merge_from(&data_from, dFrom, &data_to, dTo, &data);
511  break;
512 
513  case MODE_TO:
514  rc = c81_data_merge_to(&data_from, dFrom, &data_to, dTo, &data);
515  break;
516 
517  case MODE_BOTH:
518  rc = c81_data_merge_both(&data_from, dFrom, &data_to, dTo, &data);
519  break;
520 
521  default:
523  }
524 
525  memcpy(data.header, header, sizeof(header));
526 
527  if (rc) {
528  silent_cerr("data merge failed" << std::endl);
529  return rc;
530  }
531 
532  if (name) {
533  std::ofstream fout(name);
534  if (!fout) {
535  silent_cerr("unable to open file \"" << name << "\"" << std::endl);
537  }
538  rc = c81_data_write(fout, &data);
539 
540  } else {
541  rc = c81_data_write(std::cout, &data);
542  }
543 
544  if (rc) {
545  silent_cerr("output failed" << std::endl);
546  }
547 
548  return rc;
549 }
int NMD
Definition: aerodc81.h:160
int NAM
Definition: aerodc81.h:166
#define MBDYN_EXCEPT_ARGS
Definition: except.h:63
static int c81_data_merge_to(c81_data *data_from, doublereal wfrom, c81_data *data_to, doublereal wto, c81_data *data)
Definition: c81merge.cc:214
static int c81_data_merge_both(c81_data *data_from, doublereal wfrom, c81_data *data_to, doublereal wto, c81_data *data)
Definition: c81merge.cc:270
int NMM
Definition: aerodc81.h:165
int NML
Definition: aerodc81.h:146
static void merge_vecs(int n, doublereal *v1, doublereal w1, doublereal *v2, doublereal w2, doublereal *v)
Definition: c81merge.cc:48
doublereal * mm
Definition: aerodc81.h:167
char header[31]
Definition: aerodc81.h:144
static int c81_data_merge_from(c81_data *data_from, doublereal wfrom, c81_data *data_to, doublereal wto, c81_data *data)
Definition: c81merge.cc:107
int c81_data_read(std::istream &in, c81_data *data, const doublereal dcltol, int *ff)
Definition: c81data.cc:292
int NAL
Definition: aerodc81.h:147
doublereal * ad
Definition: aerodc81.h:163
int main(int argc, char *argv[])
Definition: c81merge.cc:279
int getopt(int argc, char *const argv[], const char *opts)
Definition: getopt.c:93
#define STRLENOF(s)
Definition: mbdyn.h:166
static void merge_data_from(int NA, int NM, int NM_to, doublereal wfrom, doublereal wto, doublereal *workv, doublereal *m_from, doublereal *m_to, doublereal *m_dst, doublereal *a_from, doublereal *a_to, doublereal *a_dst)
Definition: c81merge.cc:56
doublereal * am
Definition: aerodc81.h:168
doublereal * al
Definition: aerodc81.h:149
static const doublereal a
Definition: hfluid_.h:289
char * optarg
Definition: getopt.c:74
double doublereal
Definition: colamd.c:52
int c81_data_write(std::ostream &out, c81_data *data)
Definition: c81data.cc:610
doublereal * md
Definition: aerodc81.h:162
doublereal * ml
Definition: aerodc81.h:148
static void merge_data_to(int NA, int NM, int NM_from, doublereal wfrom, doublereal wto, doublereal *workv, doublereal *m_from, doublereal *m_to, doublereal *m_dst, doublereal *a_from, doublereal *a_to, doublereal *a_dst)
Definition: c81merge.cc:163
int NAD
Definition: aerodc81.h:161