From 8792384f5cb00cb1375da42c159c7bb8e140c497 Mon Sep 17 00:00:00 2001
From: Klaus <klaus.boehnlein@tu-dresden.de>
Date: Sat, 24 Jul 2021 16:23:05 +0200
Subject: [PATCH] reverse Load bug fixed
---
src/dune-microstructure.cc | 44 +++++++++++++++++---------------------
1 file changed, 20 insertions(+), 24 deletions(-)
diff --git a/src/dune-microstructure.cc b/src/dune-microstructure.cc
index 7fa988ff..38f9d535 100644
--- a/src/dune-microstructure.cc
+++ b/src/dune-microstructure.cc
@@ -551,7 +551,7 @@ void computeElementLoadVector( const LocalView& localView,
double energyDensity = 0.0;
for (int ii=0; ii<nCompo; ii++)
for (int jj=0; jj<nCompo; jj++)
- energyDensity += stressV[ii][jj] * forceTerm(quadPos)[ii][jj];
+ energyDensity += stressV[ii][jj] *forceTerm(quadPos)[ii][jj]; //TEST VERWENDE NEGATIVE DER LAST
size_t row = localView.tree().child(k).localIndex(i);
@@ -577,7 +577,7 @@ void computeElementLoadVector( const LocalView& localView,
double energyDensityfG = 0;
for (int ii=0; ii<nCompo; ii++)
for (int jj=0; jj<nCompo; jj++)
- energyDensityfG += stressG[ii][jj] * forceTerm(quadPos)[ii][jj];
+ energyDensityfG += stressG[ii][jj] * forceTerm(quadPos)[ii][jj]; //TEST VERWENDE NEGATIVE DER LAST
elementRhs[localPhiOffset+m] += energyDensityfG * quadPoint.weight() * integrationElement;
}
@@ -1635,8 +1635,7 @@ int main(int argc, char *argv[])
///////////////////////////////////
// Get Parameters/Data
///////////////////////////////////
-
- double gamma = parameterSet.get<double>("gamma",2.0); // ratio dimension reduction to homogenization
+ double gamma = parameterSet.get<double>("gamma",2000.0); // ratio dimension reduction to homogenization
// Material parameter laminat
// double E1 = parameterSet.get<double>("E1", 17e6); //material1
@@ -1645,8 +1644,6 @@ int main(int argc, char *argv[])
// double E2 = parameterSet.get<double>("E2", 35e6); //material2
// double nu2 = parameterSet.get<double>("nu2", 0.5);
-
-
// // Plate geometry settings
double width = parameterSet.get<double>("width", 1.0); //geometry cell, cross section
// double len = parameterSet.get<double>("len", 10.0); //length
@@ -1676,9 +1673,6 @@ int main(int argc, char *argv[])
log << "alpha: " << alpha << std::endl;
log << "gamma: " << gamma << std::endl;
-
-
-
///////////////////////////////////
// --- Choose Solver ---
@@ -1694,7 +1688,7 @@ int main(int argc, char *argv[])
///////////////////////////////////
// Generate the grid
///////////////////////////////////
- using CellGridType = YaspGrid< dim, EquidistantOffsetCoordinates< double, dim> >;
+ using CellGridType = YaspGrid<dim, EquidistantOffsetCoordinates<double, dim> >;
// Domain : [0,1)^2 x (-1/2, 1/2)
// FieldVector<double,dim> lower({0.0, 0.0, -1.0/2.0});
@@ -1704,25 +1698,22 @@ int main(int argc, char *argv[])
FieldVector<double,dim> lower({-1.0/2.0, -1.0/2.0, -1.0/2.0});
FieldVector<double,dim> upper({1.0/2.0, 1.0/2.0, 1.0/2.0});
- int nE = 10;
+ int nE = 20;
log << "Number of Elements: " << nE << std::endl;
debugLog << "Number of Elements: " << nE << std::endl;
std::array<int,dim> nElements={nE,nE,nE};
+// std::array<int,dim> nElements={100,100,1}; // TEST
// printvector(std::cout, coeffT_1 , "coeffT_1", "--" );
// printvector(std::cout, x_1 , "x_1", "--" );; //#Elements in each direction
- CellGridType grid_CE(lower,upper,nElements); //Corrector problem Domain:
+ CellGridType grid_CE(lower,upper,nElements); //Corrector Problem Domain:
using GridView = CellGridType::LeafGridView;
const GridView gridView_CE = grid_CE.leafGridView();
-
-
-
-
using ScalarRT = FieldVector< double, 1>;
using VectorRT = FieldVector< double, nCompo>;
using MatrixRT = FieldMatrix< double, nCompo, nCompo>;
@@ -1730,7 +1721,6 @@ int main(int argc, char *argv[])
using FuncScalar = std::function< ScalarRT(const Domain&) >;
using Func2Tensor = std::function< MatrixRT(const Domain&) >;
-
using VectorCT = BlockVector<FieldVector<double,1> >;
using MatrixCT = BCRSMatrix<FieldMatrix<double,1,1> >;
// using VectorCT = BlockVector<FieldVector<double,dim> >;
@@ -1770,8 +1760,6 @@ int main(int argc, char *argv[])
else
return mu2;
// return mu1;
-
-
};
@@ -1783,9 +1771,6 @@ int main(int argc, char *argv[])
};
-
-
-
auto muGridF = Dune::Functions::makeGridViewFunction(muTerm, gridView_CE);
auto muLocal = localFunction(muGridF);
auto lambdaGridF = Dune::Functions::makeGridViewFunction(lambdaTerm, gridView_CE);
@@ -1798,8 +1783,6 @@ int main(int argc, char *argv[])
-
-
/////////////////////////////////////////////////////////
// Choose a finite element space for Cell Problem
/////////////////////////////////////////////////////////
@@ -1948,6 +1931,19 @@ int main(int argc, char *argv[])
assembleCellLoad(Basis_CE, muLocal, lambdaLocal,gamma, load_alpha1 ,x3G_1);
assembleCellLoad(Basis_CE, muLocal, lambdaLocal,gamma, load_alpha2 ,x3G_2);
assembleCellLoad(Basis_CE, muLocal, lambdaLocal,gamma, load_alpha3 ,x3G_3);
+
+
+ // TEST multiply RHS with -1 :
+ for(int i=0; i< Basis_CE.size()+3 ; i++)
+ {
+ load_alpha1[i] = (-1.0)*load_alpha1[i];
+ load_alpha2[i] = (-1.0)*load_alpha2[i];
+ load_alpha3[i] = (-1.0)*load_alpha3[i];
+ }
+
+
+
+
// printmatrix(std::cout, stiffnessMatrix_CE, "StiffnessMatrix before B.C", "--");
--
GitLab