Dear S. A. Mohseni,
On 05/31/2017 01:11 PM, 'Seyed Ali Mohseni' via deal.II User Group wrote:
Dear Thomas Wick, Dear Timo Heister,
I wrote an additional postprocessor in your existing phase-field code
to allow postprocessing of strain, stress or elastic energy.
Unfortunately, it seems like the STRAIN_XX and STRAIN_XY is not
increasing in each step while the STRAIN_YY increases correctly.
This is difficult to say and heavily depends on the loading conditions.
The stress has different components, but not all
of them increase for a specific loading.
The numerical example is a simple 2d cube of length and width 1.0
consisting of 1 linear element only. The material properties are the
same from the Miehe tension example. The specimen is loaded exactly
the same way as in the tension experiment.
Of course, I am aware of the fact that the element size is totally
wrong for a phase-field simulation. I am merely trying to check the
strain in a purely elastic case before crack initiation.
I stay with my previous opinion some weeks ago: it does not make sense
to run a FE simulation only on one element !!!
Where is the problem to have 100 or 1000 elements? This simulation (see
step-3 for Laplace) also would run
only for seconds.
Is it due to the predictor-corrector nature or something which causes
this problem?
This is easy to check: just disable predictor-corrector and run with
uniform mesh refinement and then you see
whether predictor-corrector causes the problem.
But anyhow: when you work only on one element - there is no
predictor-corrector because otherwise you would have > 1 elements.
Best regards,
Thomas W.
In my solid_mechanics code written in deal.II it works correctly, the
following postprocessor implementation:
|
template<int dim>
class *Postprocessor*: public DataPostprocessor<dim>
{
public:
Postprocessor ();
void compute_derived_quantities_vector (const
std::vector<Vector<double> > &uh, const
std::vector<std::vector<Tensor<1, dim> > > &duh, const
std::vector<std::vector<Tensor<2, dim> > > &dduh, const std::vector<
Point<dim> > &normals, const std::vector<Point<dim> >
&evaluation_points, std::vector<Vector<double> > &computed_quantities)
const;
virtual std::vector<std::string> get_names () const;
virtual
std::vector<DataComponentInterpretation::DataComponentInterpretation>
get_data_component_interpretation () const;
virtual UpdateFlags get_needed_update_flags () const;
private:
void print_tensor (const Tensor<2, dim>& tensor, char* name) const;
};
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
template<int dim>
*Postprocessor*<dim>::Postprocessor ()
{
}
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
template<int dim>
void *Postprocessor*<dim>::compute_derived_quantities_vector (const
std::vector<Vector<double> > &uh, const
std::vector<std::vector<Tensor<1, dim> > > &duh, const std::vector<
std::vector<Tensor<2, dim> > > &/*dduh*/, const std::vector<Point<dim>
> &/*normals*/, const std::vector<Point<dim> > &evaluation_points,
std::vector<
Vector<double> > &computed_quantities) const
{
//TODO: Postprocessing has not been optimized for 3D yet.
// Number of quadrature points (interior)
const unsigned int n_q_points = uh.size();
// CHECK: Evaluation points
const std::vector<Point<dim> > EP = evaluation_points;
//std::cout << "\nEVALUATION POINTS\n" << std::endl;
//for (unsigned int i = 0; i < n_q_points; ++i)
//{
//for (unsigned int j = 0; j < dim; ++j)
//std::cout << " " << EP[i][j] << " ";
//std::cout << std::endl;
//}
//std::cout << std::endl;
// Constitutive matrix
SymmetricTensor<4, dim> C =
Tensors::get_elasticity_tensor<dim>(FractureMechanics<dim>::public_lambda,
FractureMechanics<dim>::public_mu);
for (unsigned int q = 0; q < n_q_points; ++q)
{
std::cout << "\n - in EVALUATION point " << q + 1 << " - \n" << std::endl;
Tensor<2, dim> grad_u;
SymmetricTensor<2, dim> eps, sigma;
double eps_ii, eps_ij, sigma_ii, sigma_ij;
// ===========================[ STRAINS ]============================
for (unsigned int i = 0; i < dim; ++i)
{
int j = i + 1;
grad_u[i] = duh[q][i];
std::cout << "DUH 0: " << duh[q][0] << std::endl;
std::cout << "DUH 1: " << duh[q][1] << std::endl;
eps = 0.5 * (grad_u + transpose(grad_u));
eps_ii = eps[i][i];
if ( j < dim )
eps_ij = eps[i][j];
//eps_ij = 2.0 * eps[i][j];
//std::cout << "STRAIN " << i << i << ": " << eps_ii << std::endl;
//std::cout << "STRAIN " << i << j << ": " << eps_ij << std::endl;
computed_quantities[q](i) = eps_ii;
computed_quantities[q](j) = eps_ij;
}
// ==========================[ STRESSES ]============================
for (unsigned int i = 0; i < dim; ++i)
{
int j = i + 1;
sigma = C * eps;
sigma_ii = sigma[i][i];
if ( j < dim )
sigma_ij = sigma[i][j];
computed_quantities[q](dim + i + 1) = sigma_ii;
computed_quantities[q](dim + j + 1) = sigma_ij;
}
// =======================[ ELASTIC ENERGY ]=========================
//double psi = 0.5 * FractureMechanics<dim>::public_lambda *
trace(eps) * trace(eps) + FractureMechanics<dim>::public_mu * eps * eps;
double psi = 0.5 * (eps[0][0] * sigma[0][0] + eps[1][1] * sigma[1][1]
+ 2.0 * eps[0][1] * sigma[0][1]);
//double psi = 0.5 * scalar_product(sigma,eps);
//std::cout << std::endl << "DISPLACEMENT GRADIENT" << std::endl;
//for (unsigned int i = 0; i < dim; ++i)
//{
//for (unsigned int j = 0; j < dim; ++j)
//std::cout << " " << std::setprecision(6) << std::fixed <<
grad_u[i][j] << " ";
//std::cout << std::endl;
//}
//std::cout << std::endl;
//print_tensor(eps, "STRAIN TENSOR");
//print_tensor(sigma, "STRESS TENSOR");
computed_quantities[q](dim * 2 + 2) = psi;
}
}
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
template<int dim>
std::vector<std::string> *Postprocessor*<dim>::get_names () const
{
std::vector<std::string> names;
std::string xyz = "xyz";
// =================[ STRAIN NAMES ]==================
// Diagonal strain components
for (int i = 0; i < dim; i++)
{
std::stringstream ss;
ss << "strain_" << xyz[i] << xyz[i];
std::string s = ss.str();
names.push_back(s);
}
// Symmetric strain components (3d may not work!)
for (int i = 0, j = 1; j < dim; j++)
{
std::stringstream ss;
ss << "strain_" << xyz[i] << xyz[j];
std::string s = ss.str();
names.push_back(s);
}
// =================[ STRESS NAMES ]==================
// Diagonal stress components
for (int i = 0; i < dim; i++)
{
std::stringstream ss;
ss << "stress_" << xyz[i] << xyz[i];
std::string s = ss.str();
names.push_back(s);
}
// Symmetric stress components (3d may not work!)
for (int i = 0, j = 1; j < dim; j++)
{
std::stringstream ss;
ss << "stress_" << xyz[i] << xyz[j];
std::string s = ss.str();
names.push_back(s);
}
// ==================[ ENERGY NAME ]==================
names.push_back("elastic_energy");
return names;
}
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
template<int dim>
std::vector<DataComponentInterpretation::DataComponentInterpretation>
*Postprocessor*<dim>::get_data_component_interpretation () const
{
std::vector<DataComponentInterpretation::DataComponentInterpretation>
interpretation(dim, DataComponentInterpretation::component_is_scalar);
interpretation.push_back(DataComponentInterpretation::component_is_scalar);
interpretation.push_back(DataComponentInterpretation::component_is_scalar);
return interpretation;
}
------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
template<int dim>
UpdateFlags *Postprocessor*<dim>::get_needed_update_flags () const
{
return update_values | update_gradients | update_quadrature_points;
}
|
Also I wonder why the elastic energy is increasing to a large value
such as 4000 or more whereby it should be around 200. But that is
another problem which only occurs in a state where the crack
propagates. Hence, when I compare the energies or strains from another
code, in the elastic regime they are almost similar, but when the
crack starts to propagate they differ.
Hope you have an idea what could be causing such observations.
Kind regards,
S. A. Mohseni
--
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see
https://groups.google.com/d/forum/dealii?hl=en
---
You received this message because you are subscribed to the Google
Groups "deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send
an email to dealii+unsubscr...@googlegroups.com
<mailto:dealii+unsubscr...@googlegroups.com>.
For more options, visit https://groups.google.com/d/optout.
--
++--------------------------------------------++
Dr. Thomas Wick
Maitre de conferences / Assistant Professor
Centre de Mathematiques Appliquees (CMAP)
Ecole Polytechnique
91128 Palaiseau cedex, France
Email: thomas.w...@cmap.polytechnique.fr
Phone: 0033 1 69 33 4579
www: http://www.cmap.polytechnique.fr/~wick/
++--------------------------------------------++
--
--
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see
https://groups.google.com/d/forum/dealii?hl=en
---
You received this message because you are subscribed to the Google Groups "deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.