Skip to content
Snippets Groups Projects
Commit 59a31b8d authored by Oliver Sander's avatar Oliver Sander Committed by sander@PCPOOL.MI.FU-BERLIN.DE
Browse files

removed more old code

[[Imported from SVN: r1716]]
parent a881d885
No related branches found
No related tags found
No related merge requests found
...@@ -785,146 +785,6 @@ assembleMatrixFD(const std::vector<Configuration>& sol, ...@@ -785,146 +785,6 @@ assembleMatrixFD(const std::vector<Configuration>& sol,
} }
template <class GridType>
void RodAssembler<GridType>::
strainDerivative(const std::vector<Configuration>& localSolution,
double pos,
Dune::FieldVector<double,1> shapeGrad[2],
Dune::FieldVector<double,1> shapeFunction[2],
Dune::array<Dune::FieldMatrix<double,2,6>, 6>& derivatives) const
{
#if 0
using namespace Dune;
assert(localSolution.size()==2);
FieldVector<double,3> r_s;
for (int i=0; i<3; i++)
r_s[i] = localSolution[0].r[i]*shapeGrad[0] + localSolution[1].r[i]*shapeGrad[1];
// Interpolate current rotation at this quadrature point
Quaternion<double> q = Quaternion<double>::interpolate(localSolution[0].q, localSolution[1].q,pos);
// Contains \parder d \parder v^i_j
array<array<FieldVector<double,3>, 3>, 3> dd_dvj;
q.getFirstDerivativesOfDirectors(dd_dvj);
Quaternion<double> q_s = Quaternion<double>::interpolateDerivative(localSolution[0].q, localSolution[1].q,
pos, 1/shapeGrad[1]);
array<Quaternion<double>,3> dq_dvj;
Quaternion<double> dq_dvij_ds[2][3];
for (int i=0; i<2; i++)
for (int j=0; j<3; j++) {
for (int m=0; m<3; m++) {
dq_dvj[j][m] = (j==m) * 0.5;
dq_dvij_ds[i][j][m] = (j==m) * 0.5 * shapeGrad[i];
}
dq_dvj[j][3] = 0;
dq_dvij_ds[i][j][3] = 0;
}
// the strain component
for (int m=0; m<3; m++) {
// the shape function
for (int i=0; i<2; i++) {
// the partial derivative direction
for (int j=0; j<3; j++) {
// parder v_m / \partial r^j_i
derivatives[m][i][j] = shapeGrad[i]*q.director(m)[j];
derivatives[m][i][j+3] = r_s *dd_dvj[m][j] * shapeFunction[i];
// \partial u_m
derivatives[m+3][i][j] = 0;
double du_dvij_i_j_m = 2* ((q.mult(dq_dvj[j])).B(m)*q_s) * shapeFunction[i][0];
Quaternion<double> tmp = dq_dvj[j];
tmp *= shapeFunction[i];
#if 1
du_dvij_i_j_m += 2 * ( q.B(m)*(q.mult(dq_dvij_ds[i][j]) + q_s.mult(tmp)));
#else
#warning Term omitted in strainDerivative
#endif
derivatives[m+3][i][j+3] = du_dvij_i_j_m;
}
}
}
#endif
}
template <class GridType>
void RodAssembler<GridType>::
rotationStrainHessian(const std::vector<Configuration>& x,
double pos,
Dune::FieldVector<double,1> shapeGrad[2],
Dune::FieldVector<double,1> shapeFunction[2],
Dune::array<Dune::Matrix<Dune::FieldMatrix<double,3,3> >, 3>& rotationDer) const
{
using namespace Dune;
assert(x.size()==2);
double eps = 1e-3;
for (int m=0; m<3; m++) {
rotationDer[m].setSize(2,2);
rotationDer[m] = 0;
}
// ///////////////////////////////////////////////////////////
// Compute gradient by finite-difference approximation
// ///////////////////////////////////////////////////////////
std::vector<Configuration> forwardSolution = x;
std::vector<Configuration> backwardSolution = x;
for (int k=0; k<2; k++) {
for (int l=0; l<3; l++) {
// infinitesimalVariation(forwardSolution[k], eps, l+3);
// infinitesimalVariation(backwardSolution[k], -eps, l+3);
forwardSolution[k].q = forwardSolution[k].q.mult(Quaternion<double>::exp((l==0)*eps,
(l==1)*eps,
(l==2)*eps));
backwardSolution[k].q = backwardSolution[k].q.mult(Quaternion<double>::exp(-(l==0)*eps,
-(l==1)*eps,
-(l==2)*eps));
Dune::array<Dune::FieldMatrix<double,2,6>, 6> forwardStrainDer;
strainDerivative(forwardSolution, pos, shapeGrad, shapeFunction, forwardStrainDer);
Dune::array<Dune::FieldMatrix<double,2,6>, 6> backwardStrainDer;
strainDerivative(backwardSolution, pos, shapeGrad, shapeFunction, backwardStrainDer);
for (int i=0; i<2; i++) {
for (int j=0; j<3; j++) {
for (int m=0; m<3; m++) {
rotationDer[m][i][k][j][l] = (forwardStrainDer[m][i][j]-backwardStrainDer[m][i][j]) / (2*eps);
}
}
}
forwardSolution = x;
backwardSolution = x;
}
}
}
template <class GridType> template <class GridType>
void RodAssembler<GridType>:: void RodAssembler<GridType>::
......
...@@ -248,23 +248,6 @@ public: ...@@ -248,23 +248,6 @@ public:
void assembleMatrixFD(const std::vector<Configuration>& sol, void assembleMatrixFD(const std::vector<Configuration>& sol,
Dune::BCRSMatrix<MatrixBlock>& matrix) const; Dune::BCRSMatrix<MatrixBlock>& matrix) const;
void strainDerivative(const std::vector<Configuration>& localSolution,
double pos,
Dune::FieldVector<double,1> shapeGrad[2],
Dune::FieldVector<double,1> shapeFunction[2],
Dune::array<Dune::FieldMatrix<double,2,6>, 6>& derivatives) const;
void rotationStrainHessian(const std::vector<Configuration>& x,
double pos,
Dune::FieldVector<double,1> shapeGrad[2],
Dune::FieldVector<double,1> shapeFunction[2],
Dune::array<Dune::Matrix<Dune::FieldMatrix<double,3,3> >, 3>& rotationDer) const;
void strainHessian(const std::vector<Configuration>& localSolution,
double pos,
Dune::array<Dune::Matrix<Dune::FieldMatrix<double,6,6> >, 3>& translationDer,
Dune::array<Dune::Matrix<Dune::FieldMatrix<double,3,3> >, 3>& rotationDer) const;
void assembleGradient(const std::vector<Configuration>& sol, void assembleGradient(const std::vector<Configuration>& sol,
Dune::BlockVector<Dune::FieldVector<double, blocksize> >& grad) const; Dune::BlockVector<Dune::FieldVector<double, blocksize> >& grad) const;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment