[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[Getfem-commits] r4556 - /trunk/getfem/src/getfem_generic_assembly.cc
From: |
Yves . Renard |
Subject: |
[Getfem-commits] r4556 - /trunk/getfem/src/getfem_generic_assembly.cc |
Date: |
Sun, 23 Mar 2014 14:57:28 -0000 |
Author: renard
Date: Sun Mar 23 15:57:27 2014
New Revision: 4556
URL: http://svn.gna.org/viewcvs/getfem?rev=4556&view=rev
Log:
some usefull nonlinear operators
Modified:
trunk/getfem/src/getfem_generic_assembly.cc
Modified: trunk/getfem/src/getfem_generic_assembly.cc
URL:
http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/getfem_generic_assembly.cc?rev=4556&r1=4555&r2=4556&view=diff
==============================================================================
--- trunk/getfem/src/getfem_generic_assembly.cc (original)
+++ trunk/getfem/src/getfem_generic_assembly.cc Sun Mar 23 15:57:27 2014
@@ -1416,6 +1416,8 @@
static ga_predef_operator_tab PREDEF_OPERATORS;
void ga_init_scalar(bgeot::multi_index &mi) { mi.resize(0); }
+ void ga_init_square_matrix(bgeot::multi_index &mi, size_type N)
+ { mi.resize(2); mi[0] = mi[1] = N; }
// Norm Operator
struct norm_operator : public ga_nonlinear_operator {
@@ -1755,7 +1757,7 @@
bool result_size(const arg_list &args, bgeot::multi_index &sizes) const {
if (args.size() != 1 || args[0]->sizes().size() != 2
|| args[0]->sizes()[0] != args[0]->sizes()[1]) return false;
- ga_init_scalar(sizes);
+ ga_init_square_matrix(sizes, args[0]->sizes()[0]);
return true;
}
@@ -1801,6 +1803,192 @@
for (size_type j = 0; j < N; ++j)
for (size_type i = 0; i < N; ++i, ++it)
*it = M(i,k)*M(l,m)*M(n,j)+M(i,m)*M(m,k)*M(l,j);
+ GMM_ASSERT1(it == result.end(), "Internal error");
+ }
+ };
+
+
+ // Right-Cauchy-Green operator (F^{T}F)
+ struct Right_Cauchy_Green_operator : public ga_nonlinear_operator {
+ bool result_size(const arg_list &args, bgeot::multi_index &sizes) const {
+ if (args.size() != 1 || args[0]->sizes().size() != 2) return false;
+ ga_init_square_matrix(sizes, args[0]->sizes()[1]);
+ return true;
+ }
+
+ // Value : F^{T}F
+ void value(const arg_list &args, base_tensor &result) const {
+ // to be verified
+ size_type m = args[0]->sizes()[0], n = args[0]->sizes()[1];
+ base_tensor::iterator it = result.begin();
+ for (size_type j = 0; j < n; ++j)
+ for (size_type i = 0; i < n; ++i, ++it) {
+ *it = scalar_type(0);
+ for (size_type k = 0; k < m; ++k)
+ *it += (*(args[0]))[i*m+k] * (*(args[0]))[j*m+k];
+ }
+ }
+
+ // Derivative : F{kj}delta{li}+F{ki}delta{lj}
+ // (comes from H -> H^{T}F + F^{T}H)
+ void derivative(const arg_list &args, size_type,
+ base_tensor &result) const { // to be verified
+ size_type m = args[0]->sizes()[0], n = args[0]->sizes()[1];
+ base_tensor::iterator it = result.begin();
+ for (size_type l = 0; l < n; ++l)
+ for (size_type k = 0; k < m; ++k)
+ for (size_type j = 0; j < n; ++j)
+ for (size_type i = 0; i < n; ++i, ++it) {
+ *it = scalar_type(0);
+ if (l == i) *it += (*(args[0]))[j*m+k];
+ if (l == j) *it += (*(args[0]))[i*m+k];
+ }
+ GMM_ASSERT1(it == result.end(), "Internal error");
+ }
+
+ // Second derivative :
+ // A{ijklop}=delta{ok}delta{li}delta{pj} + delta{ok}delta{pi}delta{lj}
+ // comes from (H,K) -> H^{T}K + K^{T}H
+ void second_derivative(const arg_list &args, size_type, size_type,
+ base_tensor &result) const { // to be verified
+ size_type m = args[0]->sizes()[0], n = args[0]->sizes()[1];
+ base_tensor::iterator it = result.begin();
+ for (size_type p = 0; p < n; ++p)
+ for (size_type o = 0; o < m; ++o)
+ for (size_type l = 0; l < n; ++l)
+ for (size_type k = 0; k < m; ++k)
+ for (size_type j = 0; j < n; ++j)
+ for (size_type i = 0; i < n; ++i, ++it) {
+ *it = scalar_type(0);
+ if (o == k) {
+ if (l == i && p == j) *it += scalar_type(1);
+ if (p == i && l == j) *it += scalar_type(1);
+ }
+ }
+ GMM_ASSERT1(it == result.end(), "Internal error");
+ }
+ };
+
+ // Left-Cauchy-Green operator (FF^{T})
+ struct Left_Cauchy_Green_operator : public ga_nonlinear_operator {
+ bool result_size(const arg_list &args, bgeot::multi_index &sizes) const {
+ if (args.size() != 1 || args[0]->sizes().size() != 2) return false;
+ ga_init_square_matrix(sizes, args[0]->sizes()[0]);
+ return true;
+ }
+
+ // Value : FF^{T}
+ void value(const arg_list &args, base_tensor &result) const {
+ // to be verified
+ size_type m = args[0]->sizes()[0], n = args[0]->sizes()[1];
+ base_tensor::iterator it = result.begin();
+ for (size_type j = 0; j < m; ++j)
+ for (size_type i = 0; i < m; ++i, ++it) {
+ *it = scalar_type(0);
+ for (size_type k = 0; k < n; ++k)
+ *it += (*(args[0]))[k*m+i] * (*(args[0]))[k*m+j];
+ }
+ }
+
+ // Derivative : F{jl}delta{ik}+F{il}delta{kj}
+ // (comes from H -> HF^{T} + FH^{T})
+ void derivative(const arg_list &args, size_type,
+ base_tensor &result) const { // to be verified
+ size_type m = args[0]->sizes()[0], n = args[0]->sizes()[1];
+ base_tensor::iterator it = result.begin();
+ for (size_type l = 0; l < n; ++l)
+ for (size_type k = 0; k < m; ++k)
+ for (size_type j = 0; j < m; ++j)
+ for (size_type i = 0; i < m; ++i, ++it) {
+ *it = scalar_type(0);
+ if (k == i) *it += (*(args[0]))[l*m+j];
+ if (k == j) *it += (*(args[0]))[l*m+i];
+ }
+ GMM_ASSERT1(it == result.end(), "Internal error");
+ }
+
+ // Second derivative :
+ // A{ijklop}=delta{ki}delta{lp}delta{oj} + delta{oi}delta{pl}delta{kj}
+ // comes from (H,K) -> HK^{T} + KH^{T}
+ void second_derivative(const arg_list &args, size_type, size_type,
+ base_tensor &result) const { // to be verified
+ size_type m = args[0]->sizes()[0], n = args[0]->sizes()[1];
+ base_tensor::iterator it = result.begin();
+ for (size_type p = 0; p < n; ++p)
+ for (size_type o = 0; o < m; ++o)
+ for (size_type l = 0; l < n; ++l)
+ for (size_type k = 0; k < m; ++k)
+ for (size_type j = 0; j < m; ++j)
+ for (size_type i = 0; i < m; ++i, ++it) {
+ *it = scalar_type(0);
+ if (p == l) {
+ if (k == i && o == j) *it += scalar_type(1);
+ if (o == i && k == j) *it += scalar_type(1);
+ }
+ }
+ GMM_ASSERT1(it == result.end(), "Internal error");
+ }
+ };
+
+
+ // Green-Lagrangian operator (F^{T}F-I)/2
+ struct Green_Lagrangian_operator : public ga_nonlinear_operator {
+ bool result_size(const arg_list &args, bgeot::multi_index &sizes) const {
+ if (args.size() != 1 || args[0]->sizes().size() != 2) return false;
+ ga_init_square_matrix(sizes, args[0]->sizes()[1]);
+ return true;
+ }
+
+ // Value : F^{T}F
+ void value(const arg_list &args, base_tensor &result) const {
+ // to be verified
+ size_type m = args[0]->sizes()[0], n = args[0]->sizes()[1];
+ base_tensor::iterator it = result.begin();
+ for (size_type j = 0; j < n; ++j)
+ for (size_type i = 0; i < n; ++i, ++it) {
+ *it = scalar_type(0);
+ for (size_type k = 0; k < m; ++k)
+ *it += (*(args[0]))[i*m+k]*(*(args[0]))[j*m+k]*scalar_type(0.5);
+ if (i == j) *it -= scalar_type(0.5);
+ }
+ }
+
+ // Derivative : (F{kj}delta{li}+F{ki}delta{lj})/2
+ // (comes from H -> (H^{T}F + F^{T}H)/2)
+ void derivative(const arg_list &args, size_type,
+ base_tensor &result) const { // to be verified
+ size_type m = args[0]->sizes()[0], n = args[0]->sizes()[1];
+ base_tensor::iterator it = result.begin();
+ for (size_type l = 0; l < n; ++l)
+ for (size_type k = 0; k < m; ++k)
+ for (size_type j = 0; j < n; ++j)
+ for (size_type i = 0; i < n; ++i, ++it) {
+ *it = scalar_type(0);
+ if (l == i) *it += (*(args[0]))[j*m+k]*scalar_type(0.5);
+ if (l == j) *it += (*(args[0]))[i*m+k]*scalar_type(0.5);
+ }
+ GMM_ASSERT1(it == result.end(), "Internal error");
+ }
+
+ // Second derivative :
+ // A{ijklop}=(delta{ok}delta{li}delta{pj} + delta{ok}delta{pi}delta{lj})/2
+ // comes from (H,K) -> (H^{T}K + K^{T}H)/2
+ void second_derivative(const arg_list &args, size_type, size_type,
+ base_tensor &result) const { // to be verified
+ size_type m = args[0]->sizes()[0], n = args[0]->sizes()[1];
+ base_tensor::iterator it = result.begin();
+ for (size_type p = 0; p < n; ++p)
+ for (size_type o = 0; o < m; ++o)
+ for (size_type l = 0; l < n; ++l)
+ for (size_type k = 0; k < m; ++k)
+ for (size_type j = 0; j < n; ++j)
+ for (size_type i = 0; i < n; ++i, ++it) {
+ *it = scalar_type(0);
+ if (o == k) {
+ if (l == i && p == j) *it += scalar_type(0.5);
+ if (p == i && l == j) *it += scalar_type(0.5);
+ }
+ }
GMM_ASSERT1(it == result.end(), "Internal error");
}
};
@@ -1936,7 +2124,12 @@
PREDEF_OPERATORS.add_method("Matrix_j1", new matrix_j1_operator());
PREDEF_OPERATORS.add_method("Matrix_j2", new matrix_j2_operator());
PREDEF_OPERATORS.add_method("Inverse", new inverse_operator());
-
+ PREDEF_OPERATORS.add_method("Right_Cauchy_Green",
+ new Right_Cauchy_Green_operator());
+ PREDEF_OPERATORS.add_method("Left_Cauchy_Green",
+ new Left_Cauchy_Green_operator());
+ PREDEF_OPERATORS.add_method("Green_Lagrangian",
+ new Green_Lagrangian_operator());
return true;
}
[Prev in Thread] |
Current Thread |
[Next in Thread] |
- [Getfem-commits] r4556 - /trunk/getfem/src/getfem_generic_assembly.cc,
Yves . Renard <=