~fdm-lab/fidlab/main

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
#include <iostream>
#include <fidlab/pde2d.hpp>

using namespace fd;
using namespace fd2;

const char *project_info =
R"(//!	Project: Nonlinear Laplace PDE with mixed Neumann-Dirichlet BCs, No. 4(6)
//! Similar to petsc2d-nl-5; this project uses 2 Neumann BCs on b1, b2.
//!
//!    y ^
//!      |
//!      |  [1] b3
//!    1 +----------+
//!      |          |
//!   [2]|          |[0]
//!   b2 |          | b1
//!      |          |
//!    0 +----------+-----> x
//!      0  [3] b4  1
//!
//! PDE:
//!	$u_{xx}+u_{yy}-uu_x=-(2+\cos(x+y))\sin(x+y)$
//!
//! BCs:
//! b1: $-\frac{\partial u}{\partial x}+u=\sin(x+y)-\cos(x+y)$
//! b2: $-\frac{\partial u}{\partial y}+u=\sin(x+y)-\cos(x+y)$
//! b3: $u=u_exact$
//! b4: $u=u_exact$
//!
//!	Exact solution:
//!	$u=\sin(x+y)$
//!
//! Boundary discretization method: does not apply (INWARD_DIFF3)
//! Nonlinear solver: PETSC
)";

int main()
{
	try
	{
		xBoundary b1(1., 0., 1.), b2(0., 1., 0.);
		yBoundary b3(1., 1., 0.), b4(0., 0., 1.);

		Mesh mesh(50, 50);
		mesh.set_boundaries({b3, b2, b4, b1});

		DirichletBC bc3(b3, "bc"), bc4(b4, "bc");
		NeumannBC bc1(b1, "-1", "1", "nbc_c1"), bc2(b2, "-1", "1", "nbc_c1");

		// Part 1: declaration of fields and functions
		Field u("u"), v("v");
		Function f("f");
		Function uex("u_exact"), bc("bc"), g("nbc_c1"), xs("exact_solution");
		FieldFunction F("nonlinear_term", 2, false, false);
		Operator Laplace("Laplace");

		// Part 2: declaration and customization of PDE
		PDE pde(&mesh);
		pde.name = "PDE2d";
		pde.boundary_discretization_method = PDE::INWARD_DIFF3;/*SYM_DIFF2*/;

		// Part 3: set PDE fields, BCs and functions
		pde.set_field(u);
		pde.set_BCs({bc1, bc2, bc3, bc4}, u);
		pde.set_functions({f, bc, g, uex, xs});
		pde.set_function(F);
	
		f.set_implementation("-(2+cos(x+y))*sin(x+y)");
		bc.set_implementation("u_exact(x,y)");
		g.set_implementation("sin(x+y)-cos(x+y)");
		uex.set_implementation("sin(x+y)");
		xs.set_implementation("u_exact(x,y)");
		
		F.set_implementation({
			/* function of variables u, v          */ "u*v",
			/* partial derivative wrt 1st variable */ "v",
			/* partial derivative wrt 2nd variable */ "u"
		});
		
		Laplace(v) = Dxx[v] + Dyy[v];

		// Part 4: set PDE equations
		pde = Laplace[u] - F(u, Dx[u]) - f;

		std::cout << "PDE";
		std::cout << (pde.linear_BCs() ? " has linear BCs" : " has nonlinear BCs") << '\n';
		std::cout << "PDE";
		std::cout << (pde.linear_eqs() ? " has linear eqs" : " has nonlinear eqs") << '\n';
		std::cout << "PDE";
		std::cout << (pde.linear() ? " is linear." : " is nonlinear.") << '\n';

		// Part 5: set methods
		pde.set_library("Petsc");

		// Part 6: output
		pde.executable = "petsc2d_nl_6";
		pde.compiler = "clang++";
		pde.project_info = project_info;
		pde.set_output("code/petsc2d_nl_6");
		pde.generate_code();
		pde.generate_makefile();
	}
	catch (const std::string& msg)
	{
		std::cout << msg << std::endl;
		exit(1);
	}
	catch (const char *msg)
	{
		std::cout << "Caught const char* exception!\n";
		std::cout << msg << std::endl;
		exit(1);
	}

	return 0;
}