~fdm-lab/fidlab/main

« back to all changes in this revision

Viewing changes to examples/gsl1d_nl_1.cpp

  • Committer: Apostol Faliagas
  • Date: 2019-02-12 21:52:42 UTC
  • Revision ID: apostol.faliagas@gmail.com-20190212215242-mua2h96aikfvog8y
Published version 0.1

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
#include <iostream>
 
2
#include "pde1d.hpp"
 
3
 
 
4
using namespace fd;
 
5
using namespace fd1;
 
6
 
 
7
const char *project_info =
 
8
R"(//!  Project: Nonlinear PDE-1d with Dirichlet BCs, No. 1
 
9
//!
 
10
//!    --+----------+-----> x
 
11
//!      0          1
 
12
//!     [b0]       [b1]
 
13
//!
 
14
//! PDE:
 
15
//!     $u_{xx}-2*u^2=0$
 
16
//!
 
17
//! BCs:
 
18
//! b0: $0$
 
19
//! b1: $1$
 
20
//!
 
21
//!     Exact solutions:
 
22
//!     $u_1(x)=\frac{3}{(2+x)^2}$
 
23
//!     $u_2(x)=?$
 
24
//!
 
25
//! Boundary discretization method: does not apply
 
26
//! Linear solver: GSL all but GMRES (it does not work).
 
27
//! Nonlinear solver: GSL all.
 
28
)";
 
29
 
 
30
int main()
 
31
{
 
32
        try
 
33
        {
 
34
                Boundary b0(0.);
 
35
                Boundary b1(1.);
 
36
 
 
37
                Mesh mesh(50);
 
38
                mesh.set_boundaries({b0, b1});
 
39
 
 
40
                DirichletBC bc0(b0, "0.75"), bc1(b1, "0.333333333333333333");
 
41
 
 
42
                // Part 1: declaration of fields and functions
 
43
                Field u("u");
 
44
                Function f("f");
 
45
                Function a("fa"), b("fb"), xs("exact_solution");
 
46
 
 
47
                // Part 2: declaration and customization of PDE
 
48
                PDE pde(&mesh);
 
49
                pde.name = "PDE1d";
 
50
                pde.boundary_discretization_method = PDE::/*SYM_DIFF2;*/INWARD_DIFF3;
 
51
 
 
52
                // Part 3: set PDE fields, BCs and functions
 
53
                pde.set_field(u);
 
54
                pde.set_BCs({bc1, bc0}, u);
 
55
                pde.set_functions({f, a, b, xs});
 
56
        
 
57
                a.set_implementation("0");
 
58
                b.set_implementation("1./3.");
 
59
                f.set_implementation("0");
 
60
                xs.set_implementation("3/((2+x)*(2+x))");
 
61
                
 
62
                // Part 4: set PDE equations
 
63
                pde = Dxx[u] - 2.*u*u;
 
64
 
 
65
                std::cout << "PDE";
 
66
                std::cout << (pde.linear_BCs() ? " has linear BCs" : " has nonlinear BCs") << '\n';
 
67
                std::cout << "PDE";
 
68
                std::cout << (pde.linear_eqs() ? " has linear eqs" : " has nonlinear eqs") << '\n';
 
69
                std::cout << "PDE";
 
70
                std::cout << (pde.linear() ? " is linear." : " is nonlinear.") << '\n';
 
71
 
 
72
                // Part 5: set methods
 
73
                pde.set_library("gsl");
 
74
                pde.set_matrix_method("LU");
 
75
                // hybrids, hybrid, dnewton, broyden, hybridsj, hybridj, newton, gnewton
 
76
                pde.set_nonlinear_solver("newton");
 
77
 
 
78
                // Part 6: output
 
79
                pde.executable = "gsl1d_nl_1";
 
80
                pde.compiler = "clang++";
 
81
                pde.project_info = project_info;
 
82
                pde.set_output("code/gsl1d_nl_1");
 
83
                pde.generate_code();
 
84
                pde.generate_makefile();
 
85
        }
 
86
        catch (const std::string& msg)
 
87
        {
 
88
                std::cout << msg << std::endl;
 
89
                exit(1);
 
90
        }
 
91
        catch (const char *msg)
 
92
        {
 
93
                std::cout << "Caught const char* exception!\n";
 
94
                std::cout << msg << std::endl;
 
95
                exit(1);
 
96
        }
 
97
 
 
98
        return 0;
 
99
}