-
Notifications
You must be signed in to change notification settings - Fork 6
/
Functions.h
145 lines (126 loc) · 4.22 KB
/
Functions.h
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
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
#ifndef FUNCTIONS_H_INCLUDED
#define FUNCTIONS_H_INCLUDED
#include <Eigen/Eigen>
#include <slepceps.h>
#include "JDMG.h"
#include "LOPGMRES.h"
typedef Eigen::Matrix<PetscScalar, -1, -1> MatrixXPS;
typedef Eigen::Matrix<PetscScalar, -1, 1> VectorXPS;
class TopOpt;
enum FUNCTION_TYPE {COMPLIANCE, VOLUME, STABILITY, FREQUENCY};
class Function_Base
{
public:
Function_Base(std::vector<PetscScalar> &values, PetscScalar min_val,
PetscScalar max_val, PetscBool objective, PetscBool calc_gradient);
virtual ~Function_Base() {}
// Compute weighted and normalized values
PetscErrorCode Compute(TopOpt *topOpt);
// Get value and gradient
PetscScalar Get_Value() {return value;}
VectorXPS &Get_Gradient() {return gradient;}
// Boolean indicating constraint or objective
PetscBool objective;
// Gradient calculation needed
PetscBool calc_gradient;
// Initialize gradient arrays
PetscErrorCode Initialize_Arrays(PetscInt nElem) {
gradients.resize(nElem, nvals); gradient.resize(nElem); return 0;
}
// Assemble all function values and gradients
static PetscErrorCode Function_Call(TopOpt *topOpt, double &f, VectorXPS &g,
VectorXPS &dgdx, MatrixXPS &dfdx);
// Print all values at the end of a run if desired
static PetscErrorCode Normalization(TopOpt *topOpt);
// Names of functions
static const char *name[];
FUNCTION_TYPE func_type;
protected:
// Number of values to compute
PetscInt nvals;
// Weight of each function value
VectorXPS weights;
// Function normalization factors
PetscScalar min_val, max_val;
// Internal function values and combined function value
VectorXPS values;
PetscScalar value;
// Internal gradients and combined gradient
MatrixXPS gradients;
VectorXPS gradient;
// Compute internal function values
virtual PetscErrorCode Function(TopOpt* topOpt) = 0;
};
class Compliance : public Function_Base
{
public:
Compliance(std::vector<PetscScalar> &values, PetscScalar min_val,
PetscScalar max_val, PetscBool objective,
PetscBool calc_gradient=PETSC_TRUE) :
Function_Base(values, min_val, max_val, objective,
calc_gradient) {func_type = COMPLIANCE;}
~Compliance() {}
protected:
PetscErrorCode Function(TopOpt *topOpt);
};
class Volume : public Function_Base
{
public:
Volume(std::vector<PetscScalar> &values, PetscScalar min_val,
PetscScalar max_val, PetscBool objective,
PetscBool calc_gradient=PETSC_TRUE) :
Function_Base(values, min_val, max_val, objective,
calc_gradient) {func_type = VOLUME;}
~Volume() {}
protected:
PetscErrorCode Function(TopOpt *topOpt);
};
class Stability : public Function_Base
{
public:
Stability(std::vector<PetscScalar> &values, PetscScalar min_val,
PetscScalar max_val, PetscBool objective,
PetscBool calc_gradient=PETSC_TRUE) :
Function_Base(values, min_val, max_val, objective,
calc_gradient), lopgmres() {
Ks = NULL; func_type = STABILITY;
}
~Stability() {MatDestroy(&Ks);}
protected:
// Stress Stiffness matrix
Mat Ks;
// Stress stiffness partial sensitivity
VectorXPS dKsdy;
// Element stress stiffness sensitivity wrt local displacement
std::vector<MatrixXPS> dKsdu;
// Adjoint vector
MatrixXPS v;
// Eigensolver
LOPGMRES lopgmres;
// Internal functions
PetscErrorCode StressFnc(TopOpt *topOpt);
MatrixXPS sigtos(VectorXPS sigma);
PetscErrorCode Function(TopOpt *topOpt);
};
class Frequency : public Function_Base
{
public:
Frequency(std::vector<PetscScalar> &values, PetscScalar min_val,
PetscScalar max_val, PetscBool objective,
PetscBool calc_gradient=PETSC_TRUE) :
Function_Base(values, min_val, max_val, objective,
calc_gradient), lopgmres() {
M = NULL; func_type = FREQUENCY;}
~Frequency() {MatDestroy(&M);}
protected:
// Mass matrix
Mat M;
// Mass matrix partial sensitivity
VectorXPS dMdy;
// Eigensolver
LOPGMRES lopgmres;
// Internal functions
PetscErrorCode DiagMassFnc(TopOpt *topOpt);
PetscErrorCode Function(TopOpt *topOpt);
};
#endif // FUNCTIONS_H_INCLUDED