« back to all changes in this revision
Viewing changes to src/pde1d.cpp
- browse files at revision 1
- compare with another revision
- download diff
- download tarball
- view history from revision 1
- Committer: Apostol Faliagas
- Date: 2019-02-12 21:52:42 UTC
- Revision ID: apostol.faliagas@gmail.com-20190212215242-mua2h96aikfvog8y
Published version 0.1
- files added:
- bin
- data
- debian
- debian/patches
- debian/source
- doc
- examples
- include
- src
added
removed
src/pde1d.cpp
1
#include <cstdlib>
2
#include <cmath>
3
#include <regex>
4
#include "pde1d.hpp"
5
#include "config_code.h"
6
#include "global_config.h"
7
8
using std::string;
9
using fd::file_exists;
10
11
// =============================================================================
12
// PDE
13
// =============================================================================
14
using fd1::PDE;
15
16
PDE::PDE(Mesh *mesh, const std::string& id) :
17
mesh_ (mesh ),
18
decl_finished_ (false),
19
defs_finished_ (false),
20
BCs_finished_ (false),
21
eqs_finished_ (false),
22
finalized_ (false),
23
library_ (GSL),
24
matrix_method_ (LU),
25
nonlinear_solver_ (DNEWTON),
26
name (id ),
27
no_jacobian_code (false),
28
boundary_discretization_method (INWARD_DIFF3)
29
{
30
if (!mesh_->finalized())
31
throw _S"PDE constructor: mesh definition is not complete.";
32
}
33
34
PDE::~PDE()
35
{
36
}
37
38
#if 0
39
void PDE::require_stage_finished(bool stage_finished, const std::string& msg)
40
{
41
if (!stage_finished)
42
}
43
44
void PDE::require_stage_not_finished(const bool& stage)
45
{
46
const bool *p = static_cast<const bool *>(&stage);
47
std::string msg;
48
if (p == &decl_finished_) msg = "";
49
if (stage_finished) throw msg;
50
}
51
#endif // 0
52
53
// -----------------------------------------------------------------------------
54
// Field manipulation
55
// -----------------------------------------------------------------------------
56
57
void PDE::set_field(const Field &fc)
58
{
59
field_.set_component(fc);
60
}
61
62
void PDE::set_fields(std::initializer_list<const Field> l)
63
{
64
for (const Field& fc : l)
65
field_.set_component(fc);
66
}
67
68
int PDE::get_component_index(const Field &fc) const
69
{
70
return field_.get_component_index(fc);
71
}
72
73
int PDE::get_component_index(int cid) const
74
{
75
return field_.get_component_index(cid);
76
}
77
78
int PDE::get_component(int index) const
79
{
80
return field_.get_component(index);
81
}
82
83
bool PDE::find_component(const Field &fc) const
84
{
85
return field_.find_component(fc);
86
}
87
88
int PDE::get_component_subindex(const Field &fc) const
89
{
90
return field_.get_component_subindex(fc);
91
}
92
93
// -----------------------------------------------------------------------------
94
// BCs
95
// -----------------------------------------------------------------------------
96
97
void PDE::set_BC(const BoundaryCondition &bc, const Field& u)
98
{
99
const int index = get_component_index(u);
100
if (index == -1) throw u.format_name() + " is not a field of " + name + ".";
101
if (find_BC(bc, u)) throw "A boundary condition was already defined "
102
"for field " + u.format_name() + _S" " +
103
"on boundary " + fd::to_string(mesh_->get_boundary_id(bc.boundary())) + ".";
104
105
switch (bc.type)
106
{
107
case BoundaryCondition::DIRICHLET:
108
set_Dirichlet_BC(bc, u);
109
break;
110
case BoundaryCondition::NEUMANN:
111
set_Neumann_BC(bc, u);
112
break;
113
case BoundaryCondition::CUSTOM:
114
set_custom_BC(bc, u);
115
break;
116
default:
117
throw "Incorrect BC specification for field " + u.format_name() +
118
" on boundary " + fd::to_string(mesh_->get_boundary_id(bc.boundary())) +
119
".";
120
}
121
}
122
123
void PDE::set_BC(const BoundaryCondition &bc)
124
{
125
THROW_IF_NOT(field_.size() == 1, "PDE must have exactly one field.");
126
set_BC(bc, field_[0]);
127
}
128
129
void PDE::set_BCs(std::initializer_list<const BoundaryCondition::Ref> l)
130
{
131
for (auto& bcref : l) set_BC(bcref.self);
132
}
133
134
void PDE::set_BCs(std::initializer_list<const BoundaryCondition::Ref> l,
135
const Field& u)
136
{
137
for (auto& bcref : l) set_BC(bcref.self, u);
138
}
139
140
bool PDE::find_BC(const BoundaryCondition& bc, const Field& u) const
141
{
142
for (BoundaryCondition *p : mesh_->BCs())
143
if (p->boundary() == bc.boundary() && p->component() == u.id) return true;
144
return false;
145
}
146
147
void PDE::set_Dirichlet_BC(const BoundaryCondition &bc, const Field& u)
148
{
149
const int index = get_component_index(u);
150
DirichletBC *dbc = new DirichletBC(static_cast<const DirichletBC&>(bc));
151
dbc->set_component(u.id);
152
mesh_->BCs().push_back(dbc);
153
const int bid = mesh_->get_boundary_id(bc.boundary());
154
switch (bid)
155
{
156
case Mesh::R: field_[index].i_end = mesh_->nx; break;
157
case Mesh::L: field_[index].i_start = 0; break;
158
default: throw _S"Internal error 0104. Please report this error.";
159
}
160
}
161
162
void PDE::set_Neumann_BC(const BoundaryCondition &bc, const Field& u)
163
{
164
const int index = get_component_index(u);
165
NeumannBC *nbc = new NeumannBC(static_cast<const NeumannBC&>(bc));
166
nbc->set_component(u.id);
167
mesh_->BCs().push_back(nbc);
168
const bool bdm2 = boundary_discretization_method == SYM_DIFF2;
169
const int bid = mesh_->get_boundary_id(bc.boundary());
170
switch (bid)
171
{
172
case Mesh::R: field_[index].i_end = bdm2 ? mesh_->nx + 1 : mesh_->nx; break;
173
case Mesh::L: field_[index].i_start = bdm2 ? -1 : 0; break;
174
default: throw _S"Internal error 0105. Please report this error.";
175
}
176
}
177
178
void PDE::set_custom_BC(const BoundaryCondition &bc, const Field& u)
179
{
180
ASSERT_MSG(bc.type == BoundaryCondition::CUSTOM, "Expected custom BC.");
181
if (boundary_discretization_method == SYM_DIFF2)
182
throw _S"Custom (user defined) boundary conditions cannot be "
183
"discretized with the traditional symmetric differences method. "
184
"Use an inward differences method.";
185
if (static_cast<const CustomBC&>(bc).expression == nullptr)
186
throw _S"Custom (user defined) boundary conditions which are undefined "
187
"cannot be selected into a PDE.";
188
189
const int index = get_component_index(u);
190
CustomBC *cbc = new CustomBC(static_cast<const CustomBC&>(bc));
191
cbc->set_component(u.id);
192
set_equation_for_BC(cbc, cbc->expression);
193
mesh_->BCs().push_back(cbc);
194
const int bid = mesh_->get_boundary_id(bc.boundary());
195
switch (bid)
196
{
197
case Mesh::R: field_[index].i_end = mesh_->nx; break;
198
case Mesh::L: field_[index].i_start = 0; break;
199
default: throw _S"Internal error 0106. Please report this error.";
200
}
201
}
202
203
// -----------------------------------------------------------------------------
204
// Set equations, functions etc
205
// -----------------------------------------------------------------------------
206
207
PDE& PDE::Ref::operator=(const Expression& expr)
208
{
209
if (eqname.empty())
210
{
211
if (eqnum < 0) throw _S"Equation numbers cannot be negative.";
212
eqname = fd::to_string(eqnum);
213
}
214
pde->set_differential_expression(eqname, expr);
215
return *pde;
216
}
217
218
PDE::EquationReplacement::EquationReplacement(int fi, int ic, double val,
219
const std::string& rv, const std::string& jv, PDE *pde) :
220
index (fi),
221
i (ic),
222
value (val),
223
rhs (rv),
224
jac (jv),
225
pde_ (pde),
226
initialized_ (false) {}
227
228
void PDE::EquationReplacement::init()
229
{
230
if (initialized_) return;
231
THROW_IF(pde_ == nullptr, "null PDE object!");
232
initialized_ = true;
233
if (!rhs.empty()) return;
234
const string fi = fd::to_string_empty_if(pde_->components() == 1, index);
235
const string eq = code_elem_lhs("e", fd::to_string(i), fi);
236
rhs = "u[" + eq + "]";
237
if (!fd::is_zero(value))
238
{
239
rhs += "-";
240
if (value > 0)
241
rhs += fd::to_string(value);
242
else
243
rhs += "("+fd::to_string(value)+")";
244
}
245
jac = "du[" + eq + "]";
246
}
247
248
void PDE::replace_equation(const Field &u, int i, const string& rhs,
249
const string& jac)
250
{
251
const int index = get_component_index(u);
252
THROW_IF(index == -1, u.format_name() + " is not a field of " + name + ".");
253
replace_equation(index, i, 0, rhs, jac);
254
}
255
256
void PDE::replace_equation(const Field &u, int i, double rhs)
257
{
258
const int index = get_component_index(u);
259
THROW_IF(index == -1, u.format_name() + " is not a field of " + name + ".");
260
replace_equation(index, i, rhs, "", "");
261
}
262
263
void PDE::replace_equation(int i, const string& rhs, const string& jac)
264
{
265
THROW_IF(field_.empty(), "There is no field in the PDE.");
266
THROW_IF(field_.size() > 1, "For multidimensional fields use "
267
"the overload that accepts a Field as its first argument.");
268
replace_equation(0, i, 0, rhs, jac);
269
}
270
271
void PDE::replace_equation(int i, double rhs)
272
{
273
THROW_IF(field_.empty(), "There is no field in the PDE.");
274
THROW_IF(field_.size() > 1, "For multidimensional fields use "
275
"the overload that accepts a Field as its first argument.");
276
replace_equation(0, i, rhs, "", "");
277
}
278
279
void PDE::replace_equation(int index, int i, double rhs_value,
280
const string& rhs,
281
const string& rhs_jac)
282
{
283
const EquationReplacement eqrepl(index, i, rhs_value, rhs, rhs_jac, this);
284
const bool present = [this, &eqrepl]() {
285
for (const auto& r : replacements_)
286
if (r.equation_is_identical_to(eqrepl)) return true;
287
return false;
288
}();
289
if (present) throw _S"The equation for field "
290
+ field_[index].format_name() + " was already replaced.";
291
replacements_.push_back(eqrepl);
292
}
293
294
void PDE::set_function(Function& f, const std::string& imp)
295
{
296
const std::string& fname = f.fname();
297
298
// Check function name
299
if (fname.empty())
300
throw _S"Function has no name.";
301
if (!fd::regex_match_funcname(fname))
302
throw fname + " is not a valid function name.";
303
304
// Check if function name is reserved
305
if (fd::find(reserved_function_names_, fname))
306
throw "The function name \'" + fname + "\' is reserved.";
307
308
// Check if function is already defined:
309
if (find_fname_in_field_functions(fname) || find_fname_in_functions(fname))
310
throw "Function " + fname + " is already selected in the PDE.";
311
312
if (!imp.empty())
313
f.set_implementation(imp);
314
315
functions_.push_back(&f);
316
}
317
318
void PDE::set_function(FieldFunction& f, const std::vector<std::string>& imp)
319
{
320
const std::string& fname = f.fname();
321
322
// Check function name
323
if (fname.empty())
324
throw _S"Field function has no name.";
325
if (!std::regex_match(fname, std::regex("[_a-zA-Z][_a-zA-Z0-9]*")))
326
throw fname + " is not a valid function name.";
327
328
// Check if function name is reserved
329
if (fd::find(reserved_function_names_, fname))
330
throw "The function name \'" + fname + "\' is reserved.";
331
332
// Check if function is already defined:
333
if (find_fname_in_field_functions(fname) || find_fname_in_functions(fname))
334
throw "Function " + fname + " is already selected in the PDE.";
335
336
if (!imp.empty())
337
f.set_implementation(imp);
338
339
field_functions_.push_back(&f);
340
}
341
342
void PDE::set_functions(std::initializer_list<Function::Ref> l)
343
{
344
for (const Function::Ref& ref : l)
345
set_function(ref.self);
346
}
347
348
void PDE::set_functions(std::initializer_list<FieldFunction::Ref> l)
349
{
350
for (const FieldFunction::Ref& ref : l)
351
set_function(ref.self);
352
}
353
354
fd1::Function& PDE::get_function(const string& func_name) const
355
{
356
for (Function *f : functions_)
357
if (f->fname() == func_name) return *f;
358
throw "Unable to find function " + func_name + " in PDE " +
359
(name.empty() ? _S"noname" : name) + ".";
360
}
361
362
fd1::FieldFunction& PDE::get_field_function(const string& func_name) const
363
{
364
for (FieldFunction *f : field_functions_)
365
if (f->fname() == func_name) return *f;
366
throw "Unable to find field function " + func_name + " in PDE " +
367
(name.empty() ? _S"noname" : name) + ".";
368
}
369
370
// -----------------------------------------------------------------------------
371
// PDE Properties
372
// -----------------------------------------------------------------------------
373
374
bool PDE::linear_BCs() const
375
{
376
for (BoundaryCondition *bc : mesh_->BCs())
377
if (bc->type == BoundaryCondition::CUSTOM &&
378
!static_cast<CustomBC *>(bc)->codeinfo.linear) return false;
379
return true;
380
}
381
382
bool PDE::linear_eqs() const
383
{
384
for (const CodeFragment& cf : equation_codes_)
385
if (!cf.linear) return false;
386
return true;
387
}
388
389
// -----------------------------------------------------------------------------
390
// Code generation utilities
391
// -----------------------------------------------------------------------------
392
393
void PDE::link(string& code) const
394
{
395
if (code.empty()) return;
396
397
bool single_component = components() == 1;
398
399
// Fix commas
400
fd::replace_all(code, "$COMMA$", single_component ? "" : ",");
401
402
// Build list if referenced fields ids
403
std::vector<int> refs;
404
std::smatch sm;
405
const string rex("\\$c\\_[1-9][0-9]*\\$");
406
auto it = code.cbegin();
407
while (std::regex_search(it, code.cend(), sm, std::regex(rex)))
408
{
409
int id;
410
try {id = std::stoi(sm.str(0).substr(3, sm.length(0) - 4));}
411
catch (std::invalid_argument& err) {
412
throw "Incorrectly formatted link: " + sm.str(0); }
413
if (!fd::find(refs, id))
414
refs.push_back(id);
415
it += sm.position(0) + sm.length(0);
416
}
417
418
for (int id : refs)
419
{
420
// Check if ref is a valid field id selected in the PDE:
421
const int cindex = get_component_index(id);
422
if (cindex == -1)
423
{
424
const string global_field_name = Field::find_global_id(id);
425
if (global_field_name.empty())
426
throw _S"Internal error 0112. An unregistered field (#" +
427
fd::to_string(id) + ") was encountered. Please report this error.";
428
throw "The field " + global_field_name + " does not belong to the PDE " +
429
_S(components() > 1 ? "system " : "") + name + ".";
430
}
431
// Fix reference:
432
const string sym = "$c_" + fd::to_string(id) + "$";
433
const string val = fd::to_string_empty_if(single_component, cindex);
434
fd::replace_all(code, sym, val);
435
}
436
437
// Check linked code
438
if (code.find('$') != string::npos) {
439
DEBUG_LOG("Generated code is possibly incorrect!");
440
} // Do not remove braces!
441
}
442
443
void PDE::link(CodeFragment& cf) const
444
{
445
link(cf.text_y);
446
link(cf.text_j);
447
link(cf.text_bx_y[0]);
448
link(cf.text_bx_y[1]);
449
link(cf.text_bx_j[0]);
450
link(cf.text_bx_j[1]);
451
}
452
453
// -----------------------------------------------------------------------------
454
// Utilities
455
// -----------------------------------------------------------------------------
456
457
void PDE::finalize()
458
{
459
finalized_ = true;
460
}
461
462
// TODO rewrite it!
463
void PDE::check() const
464
{
465
if (components() != equations())
466
throw _S"Number of equations does not match number of components.";
467
if (!finalized_) throw _S"PDE system definition is not complete.";
468
}
469
470
void PDE::check_function_implementation(bool empty, const string& name)
471
{
472
if (empty)
473
{
474
switch (fd::global_config().missing_function_impl)
475
{
476
case fd::NOT_AN_ERROR:
477
break;
478
case fd::NOTE: log_msg("NOTE: The function "+name+
479
" has no implementation. Modify the generated code file(s)"
480
" to add it manually.");
481
break;
482
case fd::WARNING: log_msg("WARNING: The function "+name+
483
" has no implementation.");
484
break;
485
case fd::EXCEPTION: throw "The function "+name+
486
" has no implementation.";
487
break;
488
}
489
}
490
}
491
492
void PDE::check_function_implementations()
493
{
494
for (const auto f: functions_)
495
check_function_implementation(f->implementation.empty(), f->fname());
496
}
497
498
void PDE::check_field_function_implementations()
499
{
500
for (const auto f: field_functions_)
501
check_function_implementation(f->implementation.empty(), f->fname());
502
}
503
504
bool PDE::boundary_discretization_method_changed() const
505
{
506
bool first_time = true;
507
int dm;
508
for (BoundaryCondition *bc : mesh_->BCs())
509
{
510
auto msg = [this, bc]() -> string
511
{
512
const int cid = bc->component();
513
const int cindex = get_component_index(cid);
514
const string field_name = cindex == -1 ? Field::find_global_id(cid) :
515
field_[cindex].format_name(true);
516
return _S" discretization method was used for the discretization "
517
+ (components() == 1 ? _S"" : "of field " + field_name + " ")
518
+ "at boundary "
519
+ fd::to_string(mesh_->get_boundary_id(bc->boundary())) + ".";
520
};
521
522
if (bc->type == BoundaryCondition::CUSTOM)
523
{
524
CustomBC *cbc = (CustomBC *) bc;
525
const int this_dm = cbc->codeinfo.boundary_discretization_method;
526
if (this_dm == -1) throw "An unknown " + msg();
527
if (this_dm == SYM_DIFF2) throw _S"Custom (user defined) boundary "
528
"conditions cannot be discretized with the traditional symmetric "
529
"differences method. Use an inward differences method.";
530
if (first_time)
531
{
532
dm = this_dm;
533
first_time = false;
534
}
535
else if (dm != this_dm)
536
{
537
throw "A different (" + fd::to_string(this_dm) + ") " + msg()
538
+ " Method previously used: " + fd::to_string(dm) + ".";
539
return true;
540
}
541
}
542
}
543
return first_time ? false : (dm != boundary_discretization_method);
544
}
545
546
bool PDE::dirichlet_neumann_BCs_on_same_boundary() const
547
{
548
for (const int bid : {Mesh::R, Mesh::L})
549
{
550
int dirichler_bcs = 0;
551
int neumann_bcs = 0;
552
for (const auto& field: field_)
553
{
554
BoundaryCondition *bc = mesh_->get_BC_for_component(field.id, bid);
555
if (bc == nullptr) continue;
556
if (bc->type == BoundaryCondition::DIRICHLET) ++dirichler_bcs;
557
else if (bc->type == BoundaryCondition::NEUMANN) ++neumann_bcs;
558
}
559
if (dirichler_bcs > 0 && neumann_bcs > 0) return true;
560
}
561
return false;
562
}
563
564
bool PDE::has_neumann_BCs() const
565
{
566
for (BoundaryCondition *bc : mesh_->BCs())
567
if (bc->type == BoundaryCondition::NEUMANN) return true;
568
return false;
569
}
570
571
bool PDE::component_has_neumann_BC_on(int comp_index, int bid) const
572
{
573
const auto bc = mesh_->get_BC_for_component(get_component(comp_index), bid);
574
THROW_IF(bc == nullptr, "no BC for component "+
575
field_[comp_index].format_name()+" on boundary "+fd::to_string(bid)+".");
576
return bc->type == BoundaryCondition::NEUMANN;
577
}
578
579
bool PDE::has_neumann_BC_on(int bid) const
580
{
581
for (int comp_index = 0; comp_index < components(); ++comp_index)
582
if (component_has_neumann_BC_on(comp_index, bid)) return true;
583
return false;
584
}
585
586
bool PDE::find_fname_in_functions(const std::string& fname) const
587
{
588
for (auto f : functions_)
589
if (fname == (f->fname())) return true;
590
return false;
591
}
592
593
bool PDE::find_fname_in_field_functions(const std::string& fname) const
594
{
595
for (auto f : field_functions_)
596
if (fname == (f->fname())) return true;
597
return false;
598
}
599
600
bool PDE::has_default_field_function_implementations() const
601
{
602
for (FieldFunction *f : field_functions_)
603
for (int n = 1; n < f->argc() + 1; ++n)
604
if ((int) f->implementation.size() > n && f->implementation[n] == "default")
605
return true;
606
return false;
607
}
608
609
bool PDE::has_exact_solution() const
610
{
611
for (auto f : functions_)
612
if (f->fname() == "exact_solution") return true;
613
return false;
614
};
615
616
// -----------------------------------------------------------------------------
617
// Select methods
618
// -----------------------------------------------------------------------------
619
620
void PDE::set_library(const std::string& library)
621
{
622
for (auto lib : libs_)
623
if (lib.second == library)
624
{ library_ = lib.first; return; }
625
throw "Unknown library: " + library;
626
}
627
628
void PDE::set_matrix_method(const std::string& matrix_method)
629
{
630
for (auto mm : matrix_methods_)
631
if (mm.second == matrix_method)
632
{ matrix_method_ = mm.first; return; }
633
throw "Unknown matrix method: " + matrix_method;
634
}
635
636
void PDE::set_nonlinear_solver(const std::string& nonlinear_solver)
637
{
638
for (auto nls : nonlinear_solvers_)
639
if (nls.second == nonlinear_solver)
640
{ nonlinear_solver_ = nls.first; return; }
641
throw "Unknown nonlinear solver: " + nonlinear_solver;
642
}
643
644
// -----------------------------------------------------------------------------
645
// Customization and options
646
// -----------------------------------------------------------------------------
647
void PDE::set_initial_estimate(const string& impl)
648
{
649
init_estimate_imp_ = impl;
650
fd::global_config().format_function_implementation(init_estimate_imp_);
651
}
652
653
bool PDE::has_initial_estimate() const
654
{ return !init_estimate_imp_.empty(); }
655
656
// -----------------------------------------------------------------------------
657
// Control output
658
// -----------------------------------------------------------------------------
659
660
std::map<const std::string, int>& PDE::get_dict_for_source(int which)
661
{
662
if (which == HEADER_FILE) return extra_headers_def_;
663
if (which == IMP_FILE) return extra_headers_imp_;
664
return extra_headers_main_;
665
}
666
667
void PDE::set_headers(const std::string& src, const std::string& pos,
668
const std::string& headers)
669
{
670
const int which = [src] {
671
if (src == "header" || src == "interface" || src == "int" || src == "h" )
672
return HEADER_FILE;
673
if (src == "implementation" || src == "imp" || src == "i") return IMP_FILE;
674
if (src == "main" || src == "m" ) return MAIN_FILE;
675
throw "Unknown output file type: " + src +".";
676
}();
677
678
const int where = [src, pos] {
679
if (pos == "top" || pos == "begin" || pos == "start" || pos == "t" )
680
return PLACEMENT_TOP;
681
if (pos == "before library" || pos == "before interface" ||
682
pos == "before header" || pos == "b" || pos == "bh" )
683
return PLACEMENT_BEFORE_LIB;
684
if (pos == "bottom" || pos == "end" ||pos == "e" ) return PLACEMENT_LAST;
685
throw "Unknown placement \"" + pos + "\" for header(s): " + src +".";
686
}();
687
688
/*std::map<const string, int>*/auto& d = get_dict_for_source(which);
689
const char delim = (headers.find(':') == string::npos) ? ';' : ':';
690
const string msg = " header specification.";
691
std::stringstream ss(headers);
692
string item;
693
694
while (std::getline(ss, item, delim))
695
{
696
if (item.empty()) throw "Incorrect" + msg;
697
if (d.find(item) != d.end()) throw "Duplicate" + msg;
698
d[item] = where;
699
}
700
}
701
702
void PDE::set_extra_includes(const std::string& s)
703
{
704
const char delim = (s.find(':') == std::string::npos) ? ';' : ':';
705
const std::string msg = " include file specification.";
706
std::stringstream ss(s);
707
std::string item;
708
while (std::getline(ss, item, delim))
709
{
710
if (item.empty()) throw "Incorrect" + msg;
711
if (fd::find(extra_includes_, item)) throw "Duplicate" + msg;
712
extra_includes_.push_back(item);
713
}
714
}
715
// -----------------------------------------------------------------------------
716
// Set output files
717
// -----------------------------------------------------------------------------
718
719
std::string PDE::get_implem_filename() const {return code_file_;}
720
std::string PDE::get_header_filename() const {return header_file_;}
721
std::string PDE::get_main_filename() const {return program_;}
722
723
void PDE::set_output(const std::string& filename, const std::string& progname)
724
{
725
auto make_filename = [](const std::string& name, bool header)->std::string {
726
static const char *ext[] = {".cpp", ".cxx", ".c++", ".cc", ".C",
727
".hpp", ".hxx", ".h++", ".hh", ".H"};
728
static const int sz[] = {4, 4, 4, 3, 2};
729
static const int n_ext = sizeof(ext)/(2*sizeof(ext[0]));
730
int w = header ? n_ext : 0;
731
if (name.rfind(ext[w + 0]) == name.size() - sz[0] ||
732
name.rfind(ext[w + 1]) == name.size() - sz[1] ||
733
name.rfind(ext[w + 2]) == name.size() - sz[2] ||
734
name.rfind(ext[w + 3]) == name.size() - sz[3] ||
735
name.rfind(ext[w + 4]) == name.size() - sz[4])
736
return name;
737
return name + ext[w];
738
};
739
// Open file makining necessary adjustments
740
auto open_file = [this, make_filename](const std::string& filename,
741
std::string& ofname, std::ofstream& of, bool hdr) {
742
ofname = make_filename(filename, hdr);
743
fd::backup_old_file(ofname);
744
of.open(ofname);
745
if (!of) throw "Could not open file " + ofname;
746
};
747
open_file(filename, code_file_, code_, false);
748
open_file(filename, header_file_, code_h_, true);
749
// Decide which file to use to output main function
750
program_ = [this, progname, make_filename]() -> std::string {
751
if (progname.empty())
752
return code_file_.substr(0, code_file_.rfind('.')) + "_main.cpp";
753
if (progname == "*") return code_file_;
754
return make_filename(progname, false);
755
}();
756
}
757
758
void PDE::switch_to_main()
759
{
760
if (program_ == code_file_) return;
761
code_.close();
762
std::string fname = get_main_filename();
763
fd::backup_old_file(fname);
764
code_.open(fname);
765
if (!code_) throw "Could not open file " + fname;
766
}
767
768
// TODO finalize PDE just before code generation
769
// TODO check that number of field components == number of eqs (in check() ?)
Loggerhead is a web-based interface for Breezy
Version: 2.0.1