?? lagrange_multiplier.cpp
字號:
#include "vs.h"
#define L_ 1.0
#define E_ 1.0
#define I_ 1.0
#define f_ 1.0
#define M_ 1.0
int main() {
{ // One parameter approximation--1. define the matrix element by element
// A. two-point Trapezoidal Integration Formula
double weight[2] = {0.5, 0.5};
Quadrature qp(weight, 0.0, L_, 2);
J d_l(L_); // the normalized length of the integration segments
// B. Define Basis Functions
H1 x(qp),
psi = x, w = x;
double lambda = 1.0;
// C. Variational Formulation
H0 b = INTEGRABLE_MATRIX("int, int, Quadrature", 3, 3, qp);
b[0][0] = (E_*I_)*d(psi)*d(psi);
b[0][2] = b[2][0] = ((H0)psi)*lambda;
b[1][2] = b[2][1] = d(w)*lambda;
b[0][1] = b[1][0] = b[1][1] = b[2][2] = 0.0;
C0 B = b | d_l; // LSH bilinear form
C0 l(3, (double*)NULL); // RHS vector
l[0] = -L_*M_; l[1] = ( ((H0)w)*f_ ) | d_l; l[2] = 0.0;
// D. Solutions--the Ritz Coefficients
C0 c = l / B;
for(int i = 0; i < 3; i++) cout << c[i] << ", "; cout << endl;
}
{ // One parameter approximation--2. use matrix concatenation
// A. two-point Trapezoidal Integration Formula
double weight[2] = {0.5, 0.5};
Quadrature qp(weight, 0.0, L_, 2);
J d_l(L_); // the normalized length of the integration segments
// B. Define Basis Functions
H1 x(qp),
psi = x, w = x;
H0 null(qp);
null = 0.0;
double lambda = 1.0;
// C. Variational Formulation
C0 B = ( ( (E_*I_*(d(psi)*d(psi)))| null | (((H0)psi)*lambda) )&
( null | null | (d(w)*lambda) )&
( (((H0)psi)*lambda) | (d(w)*lambda) | null )
) | d_l; // LSH bilinear form
C0 l = ( C0(-L_*M_) &
( (((H0)w)*f_) | d_l ) &
C0(0.0)
); // RHS vector
// D. Solutions--the Ritz Coefficients
C0 c = l / B;
for(int i = 0; i < 3; i++) cout << c[i] << ", "; cout << endl;
}
{ // One parameter approximation--3. use Cartesian basis expression
// A. two-point Trapezoidal Integration Formula
double weight[2] = {0.5, 0.5};
Quadrature qp(weight, 0.0, L_, 2);
J d_l(L_); // the normalized length of the integration segments
// B. Define Basis Functions
H1 x(qp),
psi = x, w = x;
double lambda = 1.0;
// C. Variational Formulation
C0 e(3);
C0 B = ( (E_*I_*(d(psi)*d(psi)))*(e[0]%e[0])+
0.0 *(e[0]%e[1])+
(((H0)psi)*lambda) *(e[0]%e[2])+
0.0 *(e[1]%e[0])+
0.0 *(e[1]%e[1])+
(d(w)*lambda) *(e[1]%e[2])+
(((H0)psi)*lambda) *(e[2]%e[0])+
(d(w)*lambda) *(e[2]%e[1])+
0.0 *(e[2]%e[2])
) | d_l; // LSH bilinear form
C0 l = ((-L_*M_)*e[0] + ((((H0)w)*f_) | d_l)*e[1] + 0.0*e[2]); // RHS vector
// D. Solutions--the Ritz Coefficients
C0 c = l / B;
for(int i = 0; i < 3; i++) cout << c[i] << ", "; cout << endl;
}
{ // Two parameter approximation--1. use matrix concatenation for the formulation
// A. Simpson's Integration Formula
double weight[3] = {1.0/3.0, 4.0/3.0, 1.0/3.0};
Quadrature qp(weight, 0.0, L_, 3);
J d_l(L_/2.0); // the normalized length of the integration segments
// B. Define Basis Functions
H1 x(qp),
psi = INTEGRABLE_VECTOR_OF_TANGENT_BUNDLE("int, int, Quadrature", 2/*vector size*/, 1/*spatial dim.*/, qp),
w = INTEGRABLE_VECTOR_OF_TANGENT_BUNDLE("int, int, Quadrature", 2/*vector size*/, 1/*spatial dim.*/, qp),
lambda = INTEGRABLE_VECTOR_OF_TANGENT_BUNDLE("int, int, Quadrature", 2/*vector size*/, 1/*spatial dim.*/, qp);
psi[0] = x; psi[1] = x.pow(2); // psi
w[0] = x; w[1] = x.pow(2); // w
lambda[0] = 1.0; lambda[1] = x; // lambda
// C. Weak Formulation
H0 null = INTEGRABLE_MATRIX("int, int, Quadrature", 2, 2, qp); null = 0.0;
C0 B =(( (E_*I_*(d(psi)*(~d(psi))))| null | (((H0)psi)%((H0)lambda)) )&
( null | null | (d(w)(0)%((H0)lambda)) )&
( (((H0)lambda)%((H0)psi)) | (((H0)lambda)%d(w)(0)) | null )
) | d_l;
C0 M_delta_psi(2, (double*)NULL), zero(2, (double*)NULL);
M_delta_psi[0] = -M_*L_; M_delta_psi[1] = -M_*L_*L_; zero = 0.0;
C0 l = ( M_delta_psi & ( (((H0)w)*f_) | d_l) & zero);
// D. Solutions--the Ritz Coefficients
C0 c = l / B;
for(int i = 0; i < 6; i++) cout << c[i] << ", "; cout << endl;
}
{ // Two parameter approximation--2. use referenced vector/matrix for the formulation
// A. Simpson's Integration Formula
double weight[3] = {1.0/3.0, 4.0/3.0, 1.0/3.0};
Quadrature qp(weight, 0.0, L_, 3);
J d_l(L_/2.0); // the normalized length of the integration segments
// B. Define Basis Functions
H1 x(qp),
psi = INTEGRABLE_VECTOR_OF_TANGENT_BUNDLE("int, int, Quadrature", 2/*vector size*/, 1/*spatial dim.*/, qp),
w = INTEGRABLE_VECTOR_OF_TANGENT_BUNDLE("int, int, Quadrature", 2/*vector size*/, 1/*spatial dim.*/, qp),
lambda = INTEGRABLE_VECTOR_OF_TANGENT_BUNDLE("int, int, Quadrature", 2/*vector size*/, 1/*spatial dim.*/, qp);
psi[0] = x; psi[1] = x.pow(2); // psi
w[0] = x; w[1] = x.pow(2); // w
lambda[0] = 1.0; lambda[1] = x; // lambda
// C. Weak Formulation
H0 b = INTEGRABLE_MATRIX("int, int, Quadrature", 6, 6, qp),
b00 = INTEGRABLE_MATRIX("int, int, H0&, int, int, Quadrature", 2, 2, b, 0, 0, qp), // referenced submatrix
b01 = INTEGRABLE_MATRIX("int, int, H0&, int, int, Quadrature", 2, 2, b, 0, 2, qp),
b02 = INTEGRABLE_MATRIX("int, int, H0&, int, int, Quadrature", 2, 2, b, 0, 4, qp),
b10 = INTEGRABLE_MATRIX("int, int, H0&, int, int, Quadrature", 2, 2, b, 2, 0, qp),
b11 = INTEGRABLE_MATRIX("int, int, H0&, int, int, Quadrature", 2, 2, b, 2, 2, qp),
b12 = INTEGRABLE_MATRIX("int, int, H0&, int, int, Quadrature", 2, 2, b, 2, 4, qp),
b20 = INTEGRABLE_MATRIX("int, int, H0&, int, int, Quadrature", 2, 2, b, 4, 0, qp),
b21 = INTEGRABLE_MATRIX("int, int, H0&, int, int, Quadrature", 2, 2, b, 4, 2, qp),
b22 = INTEGRABLE_MATRIX("int, int, H0&, int, int, Quadrature", 2, 2, b, 4, 4, qp);
b00 = (E_*I_)*d(psi)*(~d(psi));b01 = 0.0; b02 = ((H0)psi)%((H0)lambda);
b10 = 0.0; b11 = 0.0; b12 = d(w)(0)%((H0)lambda);
b20 = ((H0)lambda)%((H0)psi); b21 = ((H0)lambda)%d(w)(0); b22 = 0.0;
C0 B = b | d_l; // LSH bilinear form
H0 F = INTEGRABLE_VECTOR("int, Quadrature", 6, qp),
F0 = INTEGRABLE_VECTOR("int, H0&, int, const Quadrature&", 2, F, 0, qp), // referenced vector
F1 = INTEGRABLE_VECTOR("int, H0&, int, const Quadrature&", 2, F, 2, qp),
F2 = INTEGRABLE_VECTOR("int, H0&, int, const Quadrature&", 2, F, 4, qp);
F1 = (((H0)w)*f_); F0 = F2 = 0.0;
C0 l = F |d_l; // f vector size = 6
l[0] = -M_*L_; l[1] = -M_*L_*L_;
// D. Solutions--the Ritz Coefficients
C0 c = l / B;
for(int i = 0; i < 6; i++) cout << c[i] << ", "; cout << endl;
}
{ // Two parameter approximation--3. use equal-partition submatrix/subvector for the formulation
// A. Simpson's Integration Formula
double weight[3] = {1.0/3.0, 4.0/3.0, 1.0/3.0};
Quadrature qp(weight, 0.0, L_, 3);
J d_l(L_/2.0); // the normalized length of the integration segments
// B. Define Basis Functions
H1 x(qp),
psi = INTEGRABLE_VECTOR_OF_TANGENT_BUNDLE("int, int, Quadrature", 2/*vector size*/, 1/*spatial dim.*/, qp),
w = INTEGRABLE_VECTOR_OF_TANGENT_BUNDLE("int, int, Quadrature", 2/*vector size*/, 1/*spatial dim.*/, qp),
lambda = INTEGRABLE_VECTOR_OF_TANGENT_BUNDLE("int, int, Quadrature", 2/*vector size*/, 1/*spatial dim.*/, qp);
psi[0] = x; psi[1] = x.pow(2); // psi
w[0] = x; w[1] = x.pow(2); // w
lambda[0] = 1.0; lambda[1] = x; // lambda
// C. Weak Formulation
H0 b = INTEGRABLE_MATRIX("int, int, Quadrature", 6, 6, qp),
bs = INTEGRABLE_SUBMATRIX("int, int, H0&", 2, 2, b);
bs(0,0) = E_*I_*(d(psi)*(~d(psi)));bs(0,1) = 0.0; bs(0,2) = ((H0)psi)%((H0)lambda);
bs(1,0) = 0.0; bs(1,1) = 0.0; bs(1,2) = d(w)(0)%((H0)lambda);
bs(2,0) = ((H0)lambda)%((H0)psi); bs(2,1) = ((H0)lambda)%d(w)(0); bs(2,2) = 0.0;
H0 F = INTEGRABLE_VECTOR("int, Quadrature", 6, qp),
Fs = INTEGRABLE_SUBVECTOR("int, H0&", 2, F);
Fs(1) = (((H0)w)*f_); Fs(0) = Fs(2) = 0.0;
C0 B = b | d_l; // LSH bilinear form
C0 l = F | d_l; // f vector size = 6
l[0] = -M_*L_; l[1] = -M_*L_*L_;
// D. Solutions--the Ritz Coefficients
C0 c = l / B;
for(int i = 0; i < 6; i++) cout << c[i] << ", "; cout << endl;
}
{ // Two parameter approximation--4. use Cartesian basis expression for the formulation of the LHS bilinear form
// A. Simpson's Integration Formula
double weight[3] = {1.0/3.0, 4.0/3.0, 1.0/3.0};
Quadrature qp(weight, 0.0, L_, 3);
J d_l(L_/2.0); // the normalized length of the integration segments
// B. Define Basis Functions
H1 x(qp),
psi = INTEGRABLE_VECTOR_OF_TANGENT_BUNDLE("int, int, Quadrature", 2/*vector size*/, 1/*spatial dim.*/, qp),
w = INTEGRABLE_VECTOR_OF_TANGENT_BUNDLE("int, int, Quadrature", 2/*vector size*/, 1/*spatial dim.*/, qp),
lambda = INTEGRABLE_VECTOR_OF_TANGENT_BUNDLE("int, int, Quadrature", 2/*vector size*/, 1/*spatial dim.*/, qp);
psi[0] = x; psi[1] = x.pow(2); // psi
w[0] = x; w[1] = x.pow(2); // w
lambda[0] = 1.0; lambda[1] = x; // lambda
// C. Weak Formulation
C0 e(2), E(3); // Cartesian basis
C0 B =+( ((E_*I_)*(d(psi)*(~d(psi))))*((e%e)*(E[0]%E[0]))+
0.0 *((e%e)*(E[0]%E[1]))+
(((H0)psi)%((H0)lambda)) *((e%e)*(E[0]%E[2]))+
0.0 *((e%e)*(E[1]%E[0]))+
0.0 *((e%e)*(E[1]%E[1]))+
(d(w)(0)%((H0)lambda)) *((e%e)*(E[1]%E[2]))+
(((H0)lambda)%((H0)psi)) *((e%e)*(E[2]%E[0]))+
(((H0)lambda)%d(w)(0)) *((e%e)*(E[2]%E[1]))+
0.0 *((e%e)*(E[2]%E[2]))
) | d_l; // LSH bilinear form; use operator+() to convert to integrable_matrix_type
C0 M_delta_psi(2, (double*)NULL); M_delta_psi[0] = -M_*L_; M_delta_psi[1] = -M_*L_*L_;
C0 l = +( M_delta_psi *(e*E[0]) +
( (((H0)w)*f_) | d_l )*(e*E[1]) +
0.0 *(e*E[2])
);
// D. Solutions--the Ritz Coefficients
C0 c = l / B;
for(int i = 0; i < 6; i++) cout << c[i] << ", "; cout << endl;
}
return 0;
}
?? 快捷鍵說明
復制代碼
Ctrl + C
搜索代碼
Ctrl + F
全屏模式
F11
切換主題
Ctrl + Shift + D
顯示快捷鍵
?
增大字號
Ctrl + =
減小字號
Ctrl + -