main.cc 2.46 KB
Newer Older
Christoph Pflaum's avatar
Christoph Pflaum committed
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
// ------------------------------------------------------------
// main.cc
//
// ------------------------------------------------------------
#define _USE_MATH_DEFINES
#include <cmath>
#include <complex>
#include <iostream>
#include <string>
#include <vector>
#include "../program2D/source/ugblock2D.h"
#include "source/vectorVar.h"
#include "solver/iterativeSolversSystem.h"


using std::complex;
using namespace ::_COLSAMM_;

#define L2NORM(vector) sqrt(product(vector,vector).real())
#define RESTART 0


double isZero(double x)
{
	return fabs(x) < 1e-10?1.0:0.0;
}



int main(int argc, char** argv) {
	std::ofstream DATEI;

	int n = 10;

	//SectionedRectangle geo(v_x, v_y);  
        //Rectangle geo(1.0, 1.0);
	Disc geo(2.0,1.0);

	Boundary_Marker2D boundary(&geo);
	Unstructured2DGrid_Marker2D inner_points(&geo);
	inner_points.complement(&boundary);

        // boudary integral test 
        Blockgrid2D  blockBoTest(&geo,n);
        Variable2D<double> one(blockBoTest);
        Variable2D<double>  ff(blockBoTest);

        Variable2D<double>  u(blockBoTest);
        Variable2D<double>  interpol(blockBoTest);
        Variable2D<double>  g(blockBoTest);
        Variable2D<double>  p(blockBoTest);	
	
        Local_stiffness_matrix2D<double> boIntegral(blockBoTest);  

	DATEI.open("ug.vtk",std::ios :: out);
	geo.PrintVtk(&DATEI);
	DATEI.close();

	X_coordinate2d X(blockBoTest);
	Y_coordinate2d Y(blockBoTest);

/*
        Local_stiffness_matrix2D<double>  locStiff(blockBoTest);
	//	Test.Calculate(v_()*w_());
        locStiff.Calculate_HelmBoundary();        
	one = 1.0;
	ff = locStiff(one);	
	cout << "inegral Rand: " << product(ff,one) << endl;
*/

        u = X;
        g = Y -0.2;    

        VariableVector<double, Variable2D<double> > vectorU(u,u);
        VariableVector<double, Variable2D<double> > vectorB(g,g);
        VariableVector<double, Variable2D<double> > vectorI(interpol,p);  
        vectorI = vectorU + vectorB;
  
        cout << " Maximum I: " << Maximum(vectorI) << endl;
        cout << " Minimum I: " << Minimum(vectorI) << endl;

        cout << " Maximum U: " << Maximum(vectorU) << endl;
        cout << " Minimum U: " << Minimum(vectorU) << endl;

        cout << " Maximum B: " << Maximum(vectorB) << endl;
        cout << " Minimum B: " << Minimum(vectorB) << endl;
	
	
	DATEI.open("interpol.vtk");
	interpol.Print_VTK(DATEI);
	DATEI.close();	

	DATEI.open("one.vtk");
	one.Print_VTK(DATEI);
	DATEI.close();

	DATEI.open("ff.vtk");
	ff.Print_VTK(DATEI);
	DATEI.close();
}