[deal.II] Re: Problem with postprocessor

['Seyed Ali Mohseni' via deal.II User Group]

I did another trick: Increasing G_c has the same effect to achieve purely 
elastic behavior. So, I chose a high enough G_c value and it works!
I obtain the same results and the duh values grow correctly. 

Still I cannot understand fully why when cracking is initiated, there is no 
increase of the strains except the strain in the loaded direction.
Maybe you are right, I try to compare a bigger example, but I already 
compared two phase-field codes.
I have my own version with staggered scheme and it was not identical to 
your results after crack initiation.
As mentioned before, in the elastic regime everything is identical and fine.

The only thing I changed a bit is the integration from quadratic to linear, 
but I am not sure, if seeting quadrature_formula(2) alone is enough for 
linear shape functions. You also had face_quadrature_formula or this 
Lobatto integration.
Would you be so kind and explain what lines I have to change in order to 
obtain linear shape functions and not quadratic like set in the original 
version?
How many changes and the code would help a lot.
Just to be sure I did it correct.

Thank you.

Best,
Seyed Ali 

-- 
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.

[deal.II] Re: Problem with postprocessor

['Seyed Ali Mohseni' via deal.II User Group]

Dear Thomas,

I already run the example with uniform mesh, hence global refinement with 0 
refinement cycles.
The problem is, a specimen with 100 or 1000 elements is difficult to check 
due to the terminal output being flooded.
I could just output the result for one element, but then wheres the 
difference. 
Don't you ever do some simple FE patch examples?

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Problem with postprocessor

['Seyed Ali Mohseni' via deal.II User Group]

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.
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.   

Is it due to the predictor-corrector nature or something which causes this 
problem?

In my solid_mechanics code written in deal.II it works correctly, the 
following postprocessor implementation:

template
class *Postprocessor*: public DataPostprocessor
{

public:

Postprocessor ();

void compute_derived_quantities_vector (const std::vector 
, const std::vector > > , const 
std::vector > > , const std::vector<
Point > , const std::vector _points, 
std::vector _quantities) const;

virtual std::vector get_names () const;

virtual 
std::vector 
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
*Postprocessor*::Postprocessor ()
{
}

--

template
void *Postprocessor*::compute_derived_quantities_vector (const 
std::vector , const std::vector > > , const std::vector<
std::vector > > &/*dduh*/, const std::vector 
&/*normals*/, const std::vector _points, 
std::vector<
Vector > _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 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(FractureMechanics::public_lambda, 
FractureMechanics::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::public_lambda * trace(eps) * 
trace(eps) + FractureMechanics::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;
}
}

[deal.II] Re: Crack propagation

['Seyed Ali Mohseni' via deal.II User Group]

Dear Denis Davydov,

I am trying to implement the pack_values() and unpack_values() functions 
for a FullMatrix, but the function itself in 
quadrature_point_data.h is defined using a vector of double values only. 
Nevertheless, my initial approach would be the following:

struct MyQData: public TransferableQuadraturePointData
{
double psi;
FullMatrix B;

unsigned int number_of_values () const
{
return 2;
}

virtual void pack_values (std::vector ) const
{
Assert(scalars.size()==1, ExcInternalError());
scalars[0] = psi;
}

virtual void pack_values (std::vector ) const
{
Assert(scalars.size()==1, ExcInternalError());
matrices[0] = B;
}

virtual void unpack_values (const std::vector )
{
Assert(scalars.size() ==1, ExcInternalError());
psi = scalars[0];
}

virtual void unpack_values (const std::vector )
{
Assert(matrices.size() ==1, ExcInternalError());
B = matrices[0];
}
};

Is the above correct? Will it work during refinement? I checked it without 
refinement and it seems to give me correct values at each Gauss point after 
data storage.


Note that this may not be what you want in the context of crack propagation 
> because what we do to transfer is assume that the quadrature data is 
> continuous across the element (thus the name 
> ContinuousQuadratureDataTransfer).
>

Yes, but what choice do I have? Is there another computationally efficient 
approach? 
This means probably that the PointHistory function data storage based on 
step-18 has the same drawback for discontinuous problems, am I right?  

I would be gladly willing to implement the 
DiscontinuousQuadratureDataTransfer function, if it is somehow possible to 
accomplish in deal.II. I am pretty sure it will be beneficial for others 
also in many ways. ;) 

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Crack propagation

['Seyed Ali Mohseni' via deal.II User Group]

Dear Thomas,

Sure, my bad. I am sorry for the inconvenience it may have caused you and 
Timo. Unfortunately, I was just paying attention to the input text file and 
oversaw the changes I actually made. Won't happen again. ;)

Regarding the speedup, I just figured out how to make use of the 
CellDataStorage function which has been recently implemented within 
deal.II. First of all, I want to thank Denis Davydov and Jean-Paul Pelteret 
for implementing it so nice and stable. :) Without their recent work, I 
probably would have to find a way to write it myself somehow. Regarding the 
CellDataStorage approach: Is it possible to store a FullMatrix 
within the MyQData struct function? If yes, do I have to write separate 
functions for pack_values and unpack_values functions?

Kind regards,
Seyed Ali

-- 
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.

[deal.II] Re: Crack propagation

['Seyed Ali Mohseni' via deal.II User Group]

Hi,

I just figured out how to make use of the CellDataStorage function which 
has been recently implemented within deal.II. First of all, I want to thank 
Denis Davydov and Jean-Paul Pelteret for implementing it so nice and 
stable. :) Without your recent work, I probably would have to find a way to 
write it myself somehow. Regarding the CellDataStorage approach: Is it 
possible to store a FullMatrix within the MyQData struct function? 
If yes, do I have to write separate functions for pack_values and 
unpack_values functions?

Kind regards,
Seyed Ali

-- 
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.

[deal.II] Re: Crack propagation

['Seyed Ali Mohseni' via deal.II User Group]

Dear Thomas,

I found the problem. Your version works exactly like discussed in your 
paper, even faster with my Intel i7, namely around 1400 s. The 
timeconsuming part is the computation of my B-operator.


The B-operator computation is implemented inside* assemble_system (bool 
residual_only):*

for (; cell != endc; ++cell)
if ( cell->is_locally_owned() )
{
fe_values.reinit(cell);

local_matrix = 0;
local_rhs = 0;

int cell_index = cell->active_cell_index();

...

// COMPUTE: B-Operator

for (unsigned int q = 0; q < n_q_points; ++q) // loop over 
integration/quadrature points (GPs)
{
 b_operators[cell_index][q] = Tensors::*get_b_operator*(fe_values, 
dofs_per_cell, q);
}

*B* = b_operators;

// Old Newton iteration values
fe_values.get_function_values(rel_solution, old_solution_values);
fe_values.get_function_gradients(rel_solution, old_solution_grads);

...
  
}


The variable *B *is a history variable I created to store the B-operator 
history and has been initialized in the main class function according to:

std::vector>B; // temporary cell storage 
solution (inefficient!)


The function *get_b_operator()* which is called within the Gauss points 
loop is implemented in Tensors:

template
FullMatrix get_b_operator (const FEValues _values, const 
unsigned int dofs_per_cell, const unsigned int q)
{
FullMatrix tmp(dim, GeometryInfo::vertices_per_cell);

// Remark: For vector-valued problems each column is a value for each value 
of the solution variable, hence here 3 DoFs per node (dofs_per_cell = 12)!
for (unsigned int i = 0; i < dofs_per_cell; i += 3)
{
const unsigned int index = i / 3;

// COMPUTE: B-operator (Remark: This version has to be extended for 3D!)
tmp[0][index] = fe_values.shape_grad_component(i, q, 0)[0];
tmp[1][index] = fe_values.shape_grad_component(i, q, 0)[1];
}

return tmp;
}


I assume (as I already wrote within my comments inside the above code) that 
this kind of cell storage for matrices turns out to be quite slow and 
inefficient. Jean-Paul offered me once to use the new CellDataStorage 
function, if I remember correctly. Isn't there a nice and elegant way to 
store data in each Gauss point and cell without loss of computational time? 
Or am I able to pass the B-operator to another function somehow after I 
have solved the system? Something similar to your predictor-corrector 
scheme. It seems you somehow pass data more efficiently from one function 
to another function even after you solved the system. Do you also pass cell 
data from a function before solving to a function after solving?

Best,
Seyed Ali



-- 
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.

[deal.II] Re: Crack propagation

['Seyed Ali Mohseni' via deal.II User Group]

Hi,
 

> Each Core i7 CPU has maximum 4 physical cores. So, in MPI mode, solving 
> the simulation using 4 core is faster than 8 cores in your computer.
>

This is not true. I checked it, it is slower.

Best,
Seyed Ali 

-- 
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.

[deal.II] Re: Crack propagation

['Seyed Ali Mohseni' via deal.II User Group]

Dear Thomas,

I was wondering, why my computation is rather slow compared to your results 
within your paper. For instance regarding the shear test you mention that 
for global refinement the total time amounts to 5036 s. In my case each 
step takes around 200 s which would result in a total time of 30200 s (for 
151 steps). What could be the reason?

Attached you find the terminal output from your approach.



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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Crack propagation

['Seyed Ali Mohseni' via deal.II User Group]

Dear Thomas,

Thank you. I checked the values and they are similar.

Within your predictor-corrector approach, I assume that you somehow 
exchange data between each step in order to be able to interpolate the old 
solution. Am I able to use the same functions to store energy data in each 
Gauss point within each element? Or am I forced to implement and use the 
quadrature point history approach from step-18? Due to the storage of all 
the values in each element, the computation becomes inefficient. It would 
be nice to be able to access data anywhere and anytime. Actually, if you 
are using Miehe's approach I am confident that you will have to store the 
energy somewhere in order to find the maximum value of the energy. Is that 
true?

Kind regards,
Seyed Ali 

-- 
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.

[deal.II] Re: Crack propagation

['Seyed Ali Mohseni' via deal.II User Group]

Dear Thomas,

Exactly, I did that already. It seems it is identical, but I would like to 
see, if your energy variable psi_e is identical to mine. This is quite 
cumbersome for 500 elements or more to check. I mean the computation in the 
first steps is anyway independent from the phase-field value due to being 
in the elastic regime and the body has not been fractured yet. So, I assume 
it should not be a problem to just check the energy within the first steps. 
Of course you are right, the moment we reach the fracture toughness G_c 
cracks will start to grow and we have to separate between crack energy and 
bulk energy.

I cite your paper ;)


The results illustrate the expected behavior: as long as the crack does not 
grow (up to t = 0.0095 s), there is only an increase in bulk energy. Once 
the crack starts growing, bulk energy is dissipated into crack energy.


By the way Thomas, if you don't remember me, we already met after your GACM 
2015 conference presentation and within the GAMM phase-field workshop :) 
Maybe you remember: Back then you already mentioned that h < epsilon to 
resolve the transition of the phase-field variable. Nevertheless, this is 
true for the cracked region, but if you have only elastic case the mesh 
size will not matter. This is just to allow a correct and accurate 
resolution and regularization of the phase-field crack region. 

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Crack propagation

['Seyed Ali Mohseni' via deal.II User Group]

Dear Thomas,

I would like to check the elastic energy you compute and compare it to my 
results. For a single element the computation in an elastic regime should 
be independent from the phase-field approach in primary steps.

BR,
Seyed Ali 

-- 
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.

[deal.II] Re: Crack propagation

['Seyed Ali Mohseni' via deal.II User Group]

Dear Thomas,

I tried what you suggested and it works. Thank you :)
There is a slight problem though when I try to compute the same example 
with the same boundary conditions for the unit_square_4.inp file. I am 
trying to run a one element test using your phase-field approach. It works, 
but there is no residual value and next to it I obtain a -nan below 
Reduction. Is it possible to run a one element test with the same settings 
or do I have to modify something (such as rhs computations or boundary 
value functions)?

Kind regards,
Seyed Ali

-- 
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.

[deal.II] Crack propagation

['Seyed Ali Mohseni' via deal.II User Group]

Dear Timo Heister and Thomas Wick,

I am trying to run your phase-field crack propagation example, namely the 
single edge notched tension test, with less than 9 Gauss points. 
Hence, the fe variable is initialized using FE_Q(degree) where degree 
cannot be 0. This means degree should always be >1 here in order to make 
things work, but how can one compute the examples with 1 or 4 Gauss points? 
Is the code not suited for what I am trying to achieve or where are the 
parts I have to change in order to make it work for less integration 
points. To be honest, I don't even see the point of having 9 integration 
points for a 4-noded quadrilateral or am I overseeing something?

Thank you in advance.

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Interaction with the boundary description for a parallel distributed triangulation

['Seyed Ali Mohseni' via deal.II User Group]

Dear Jean-Paul,
 

> Yes, this looks like a copy-paste error in the documentation. Would you 
> like to submit a patch to correct it? 
>

Sure. How and where can I submit a patch?  

When you refine cells at a boundary, the child cell's faces will inherit 
> the boundary and manifold ID's of the parent cells regardless of whether 
> the mesh was read in or created using the GridGenerator functions. 
>

But according to the documentation this is only true for sequential 
triangulations and not for distributed triangulations. It is more safe to 
update boundary IDs while the 
distributed::parallel::Triangulation::execute_coarsening_and_refinement 
function is running.

I cite: 

This is not necessary for sequential triangulations because, there, these 
flags are inherited from mother to child cell and remain with a cell even if 
it is refined and the children are later coarsened again, but this does not 
hold for distributed triangulations.

 See: 
https://www.dealii.org/developer/doxygen/deal.II/classparallel_1_1distributed_1_1Triangulation.html

There is a slight problem using the approach written there. I am trying to 
compile the following:

triangulation.signals.post_refinement.connect(std_cxx11::bind(&
SolidMechanics::set_boundary_ids, std_cxx11::cref(*this), std_cxx11::
ref(triangulation)));

Then I receive the following compilation error:

/usr/include/c++/5/functional(1426): error: static assertion failed with 
"Wrong number of arguments for pointer-to-member"
static_assert(_Varargs::value
^
  detected during:
instantiation of class "std::_Bind_check_arity<_Tp _Class::*, 
_BoundArgs...> [with _Tp=void (), _Class=SolidMechanics<2>, 
_BoundArgs=, 
std::reference_wrapper>>]" at line 1440
instantiation of class "std::_Bind_helper<_SocketLike, _Func, 
_BoundArgs...> [with _SocketLike=false, _Func=void 
(SolidMechanics<2>::*)(), _BoundArgs=, 
std::reference_wrapper>>]" at line 787 of 
"/home/seyedali/fe_models/deal.II/solid_mechanics/solid_mechanics.cc"
instantiation of "void SolidMechanics::create_geometry() 
[with dim=2]" at line 544 of 
"/home/seyedali/fe_models/deal.II/solid_mechanics/solid_mechanics.cc"
instantiation of "void 
SolidMechanics::do_initial_timestep() [with dim=2]" at line 1440 of 
"/home/seyedali/fe_models/deal.II/solid_mechanics/solid_mechanics.cc"
instantiation of "void SolidMechanics::run() [with dim=2]" 
at line 1472 of 
"/home/seyedali/fe_models/deal.II/solid_mechanics/solid_mechanics.cc"

/home/seyedali/programming/c++/libraries/boost-1.62.0/include/boost/function/function_template.hpp(159):
 
error: no instance of overloaded function "std::_Bind<_Functor 
(_Bound_args...)>::operator() [with _Functor=std::_Mem_fn::*)()>, _Bound_args=, 
std::reference_wrapper>>]" matches the argument list
object type is: std::_Bind (std::reference_wrapper>, 
std::reference_wrapper>)>
BOOST_FUNCTION_RETURN((*f)(BOOST_FUNCTION_ARGS));
^
/usr/include/c++/5/functional(1171): note: this candidate was rejected 
because at least one template argument could not be deduced
  operator()(_Args&&... __args) const volatile
  ^
/usr/include/c++/5/functional(1157): note: this candidate was rejected 
because at least one template argument could not be deduced
  operator()(_Args&&... __args) volatile
  ^
/usr/include/c++/5/functional(1143): note: this candidate was rejected 
because at least one template argument could not be deduced
  operator()(_Args&&... __args) const
  ^
/usr/include/c++/5/functional(1129): note: this candidate was rejected 
because at least one template argument could not be deduced
  operator()(_Args&&... __args)
  ^
  detected during:
instantiation of "void 
boost::detail::function::void_function_obj_invoker0::invoke(boost::detail::function::function_buffer &) [with 
FunctionObj=std::_Bind 
(std::reference_wrapper>, 
std::reference_wrapper>)>, R=void]" at line 936
instantiation of "void boost::function0::assign_to(Functor) 
[with R=void, Functor=std::_Bind (std::reference_wrapper>, 
std::reference_wrapper>)>]" at line 727
instantiation of "boost::function0::function0(Functor, 
boost::enable_if_c<, int>::type) [with R=void, 
Functor=std::_Bind 
(std::reference_wrapper>, 
std::reference_wrapper>)>]" at line 1073
instantiation of "boost::function::function(Functor, 
boost::enable_if_c<, int>::type) 

[deal.II] Interaction with the boundary description for a parallel distributed triangulation

['Seyed Ali Mohseni' via deal.II User Group]

Hi,

I assume there is a mistake within the documentation of 
parallel::distributed::Triangulation class. Below the 
"detailed description" part within the second box there is written 
. Doesn't it have to be ?

With regard to the same topic: I would like to know, if it is possible to 
read from a .msh input file and refine afterwards for a parallel 
distributed mesh. Hence, using the procedure from the documentation by 
means of the set_boundary_ids() function there should be a way to reuse the 
predefined boundary ID structure from the input file, am I right? How can 
we reuse the initial boundary ID setting efficiently without much effort?  

In my case mesh refinement works for 1 core, but for several cores it fails 
with the following error:

Assertion failed in file 
../src/mpid/ch3/channels/nemesis/src/ch3_progress.c at line 560: 
!vc_ch->recv_active
internal ABORT - process 0
Assertion failed in file 
../src/mpid/ch3/channels/nemesis/src/ch3_progress.c at line 560: 
!vc_ch->recv_active
internal ABORT - process 4

I assume this problem rises from the boundary description since when I 
refine adjacent cells to other nodes not directly related to the boundary, 
it also works in parallel mode.

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Access specific element within a distributed triangulation

['Seyed Ali Mohseni' via deal.II User Group]

Dear Bruno,

Thanks. That's even better and even in just a single line of code. ;)

MPI_Bcast(, 1, MPI_DOUBLE, max_rank, mpi_com);

Best,
Seyed Ali

-- 
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.

[deal.II] Re: Access specific element within a distributed triangulation

['Seyed Ali Mohseni' via deal.II User Group]

Dear Bruno,

Thanks. That's even better and even in just a single line of code. ;)

MPI_Bcast(, 1, MPI_DOUBLE, max_rank, mpi_com);



Best,
Seyed Ali

-- 
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.

[deal.II] Re: Access specific element within a distributed triangulation

['Seyed Ali Mohseni' via deal.II User Group]

Dear all,

I was able to solve the MPI variable problem. 

For the users who may have the same issue once:

if ( this_mpi_process == max_rank )
{
 for (unsigned int i = 0; i < n_mpi_processes; ++i)
 {
  if ( i != max_rank )
   MPI_Send(, 1, MPI_DOUBLE, i, 0, mpi_com);
 }
}
else
{
 MPI_Recv(, 1, MPI_DOUBLE, max_rank, 0, mpi_com, 
MPI_STATUS_IGNORE);
 printf("Process %d received number from process %d\n", 
this_mpi_process, max_rank);
}

A double variable called position (just as an example here, you have to 
pass x,y,z, so 3 variables) can be defined earlier which has then to be 
filled on the rank equal to max_rank. 

Finally, I can proceed. The code was tested for ranks 1 to 16 and works 
perfectly. Every time the maximum value and position remains correct and 
unchanged, but the maximum rank (max_rank) changes dependent on the core 
numbers.

Thank you a lot. 
The community here is extraordinary. ;)

Best 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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Access specific element within a distributed triangulation

['Seyed Ali Mohseni' via deal.II User Group]

Dear all,

I was able to solve the MPI variable problem. 

For the users who may have the same issue once:

if ( this_mpi_process == max_rank )
{
 for (unsigned int i = 0; i < n_mpi_processes; ++i)
 {
  if ( i != max_rank )
   MPI_Send(, 1, MPI_DOUBLE, i, 0, mpi_com);
 }
}
else
{
 MPI_Recv(, 1, MPI_DOUBLE, max_rank, 0, mpi_com, 
MPI_STATUS_IGNORE);
 printf("Process %d received number from process %d\n", 
this_mpi_process, max_rank);
}

A double variable called position (just as an example here, you have to 
pass x,y,z, so 3 variables) has been defined earlier which was filled on 
the rank equal to max_rank. 

Finally, I can proceed. The code was tested for ranks 1 to 16 and works 
perfectly. Every time the maximum value and position remains correct and 
unchanged, but the maximum rank (max_rank) changes dependent on the core 
numbers.

Thank you a lot. 
The community here is extraordinary. ;)

Best 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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Access specific element within a distributed triangulation

['Seyed Ali Mohseni' via deal.II User Group]

@Prof. Bangerth: I will try to learn debugging parallel code and also to 
follow your words of advice. Please do not apologize since I don't think 
you ever became rude. I mean, after 1 emails I would definitely not be 
as calm as you. Thank you.  :) 

@Daniel: Thumbs up. From your words I realize, there is still a lot to 
learn about deal.II yet. To be honest, I always try to solve it in a simple 
way than I extend it to the general case such as 3D, DG elements, etc. 
Thank you for all the help.

Extra question: Can we store variables or copy them to all processors? 
Since I am filling a variable in a locally owned cell, currently on rank 3. 
Then my other ranks, especially the root rank 0 has no clue of the values 
which are set. Hence, there has to be a possibility to copy data to other 
processors or let them know about the content?

Kind regards,
Seyed Ali
 

-- 
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.

[deal.II] Re: Access specific element within a distributed triangulation

['Seyed Ali Mohseni' via deal.II User Group]

Dear all,

With MPI_Allreduce it works like a charm. Thank you very much everyone.

Especially Prof. Bangerth and Daniel.

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Access specific element within a distributed triangulation

['Seyed Ali Mohseni' via deal.II User Group]

I think I figured it out ;) 
After thinking about your suggestion, Prof. Bangerth.
It came to my mind, that the result of my max_rank variable is stored on a 
specific processor due to MPI_Reduce.
That means I cannot access these data parts without being currently at the 
correspinding processor rank.
For instance I store my results on rank 0, if I want to check the max_rank 
data on processor rank 3, it is empty.
That is why I cannot get any output, because when rank 3 owns everything 
output works fine.

Now my question how can I store the data on all processors?
Or how am I able to at least store my max_rank variable on all processors?

BR,
Seyed Ali 

-- 
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.

[deal.II] Re: Access specific element within a distributed triangulation

['Seyed Ali Mohseni' via deal.II User Group]

Dear Daniel,

Thanks a lot. I understand now.
Example 3 is exactly what I need, I try the MPI_MAXLOC approach.

But is it necessary to #include  ?
Since deal.II won't recognize functions such as MPI_Comm_rank 
(MPI_COMM_WORLD, 
) for instance or MPI_Reduce.
Is it implemented in Utilities::MPI::something ? 

Kind regards,
Seyed Ali

-- 
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.

[deal.II] Re: Access specific element within a distributed triangulation

['Seyed Ali Mohseni' via deal.II User Group]

Dear Daniel,

Thank you. Does your approach also apply to MPICH? Since I am using MPICH, 
I wonder, if this command can be used there also.
This struct variable, is there an example within deal.II steps? Or do you 
have an example in C++ I can learn from. 

To be honest I am not yet very familiar with MPI itself, barely 
understanding deal.II parallelism ;)

I have also another idea: How about I find the maximum value regardless of 
the position and processor core. Then just loop again over all cells and 
compare the maximum value with the existing value within each cell. Would 
this be a very inefficient approach or is looping over cells usually 
computationally cheap in deal.II.

Kind regards,
Seyed Ali 

-- 
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.

[deal.II] Re: Access specific element within a distributed triangulation

['Seyed Ali Mohseni' via deal.II User Group]

Dear Prof. Bangerth and Jean-Paul,

What I understood so far is, that the Utilities::MPI::max() function 
computes the maximum between all processors for each row of my vector. This 
means if I have a vector of 20 entries, each entry has 8 values distributed 
on 8 cores for instance. As a result, I get 20 maximum values again, but 
this time they are maximal and processor independent. 

Hence, my question now is how am I able to check where the maximum value 
comes from? Is it possible to output the processor number of the maximum 
value for each row?
My aim is to find the node where the maximum value exists. This is somehow 
cumbersome in parallel mode since I have to store information about the 
current CPU core, cell ID and node number. Then somehow use my previously 
shown distance function in C++ and check the position of the maximum value 
with regard to the geometry structure I stored.

Is there a more elegant way to solve my problem than what I suggested in 
deal.II?

Thank you for your kind assistance so far :)
And hopefully you overlook my style of approaching things in this rather 
silly and questionable way.

Kind regards,
Seyed Ali

-- 
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.

[deal.II] Re: Access specific element within a distributed triangulation

['Seyed Ali Mohseni' via deal.II User Group]

Thank you everyone. I understand your suggestions. From what I learned now 
in deal.II is, that each cell is owned by a specific processor which means 
I cannot access the information within a locally owned cell. That's why I 
have to loop over all cells and work with distributed data without using 
locally owned entities. 

Maybe I describe my aim more properly. I like to find the maximum value 
within a MPI vector such as our locally_relevant_solution in step-40. 
Hence, I have to loop over all cells and all nodes to find where the 
geometrical position and value of the displacement lies. How is this 
possible within a parallely distributed environment?

My approach using local_dof_indices is somehow impossible. It only works 
for 1 core according to the below code:

std::vector displacement_norms, loaded_nodes, loaded_elements;
double max_displacement;
int max_loaded_node, max_loaded_element;

typename DoFHandler::active_cell_iterator cell = 
dof_handler.begin_active(), endc = dof_handler.end();
std::vector local_dof_indices(fe.dofs_per_cell);
unsigned int cell_number = 0;

for (; cell != endc; ++cell, ++cell_number)
{
 pcout << "CELL ID: " << cell_number << std::endl; // this works 
surprisingly.

 cell->get_dof_indices(local_dof_indices);// Output or accessing 
this won't work I assume.

 for (unsigned int i = 0, j = 1; i < fe.dofs_per_cell; i += 2, j += 2)
 {
  double displacement_x = locally_relevant_solution(local_dof_indices[i]);
  double displacement_y = locally_relevant_solution(local_dof_indices[j]);

  double displacement_norm = sqrt(displacement_x * displacement_x + 
displacement_y * displacement_y);

  loaded_nodes.push_back(cell->vertex_index(v));
  loaded_elements.push_back(cell_number);
  displacement_norms.push_back(displacement_norm);
 }
}

max_displacement = *max_element(displacement_norms.begin(), 
displacement_norms.end());

double max_index = distance(displacement_norms.begin(), max_element(
displacement_norms.begin(), displacement_norms.end()));

// Node and element number of maximal displacements
max_loaded_node = loaded_nodes[max_index];
max_loaded_element = loaded_elements[max_index];


Hopefully someone has an idea to help me out. Thank you :)


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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Access specific element within a distributed triangulation

['Seyed Ali Mohseni' via deal.II User Group]

Thank you everyone. I understand your suggestions. From what I learned now 
in deal.II is, that each cell is owned by a specific processor which means 
I cannot access the information within a locally owned cell. That's why I 
have to loop over all cells and work with distributed data without using 
locally owned entities. 

Maybe I describe my aim more properly. I like to find the maximum value 
within a MPI vector such as our locally_relevant_solution in step-40. 
Hence, I have to loop over all cells and all nodes to find where the 
geometrical position and value of the displacement lies. How is this 
possible within a parallely distributed environment?

My approach using local_dof_indices is somehow impossible:

for (; cell != endc; ++cell, ++cell_number)
{
 pcout << "CELL ID: " << cell_number << std::endl; // this works 
surprisingly.

 cell->get_dof_indices(local_dof_indices);// Output or accessing 
this won't work I assume.

 for (unsigned int i = 0, j = 1; i < fe.dofs_per_cell; i += 2, j += 2)
 {
  double displacement_x = locally_relevant_solution(local_dof_indices[i]);
  double displacement_y = locally_relevant_solution(local_dof_indices[j]);

double config_forces_norm = sqrt(config_forces_x * config_forces_x + 
config_forces_y * config_forces_y);
 }
}


-- 
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.

[deal.II] Re: Access specific element within a distributed triangulation

['Seyed Ali Mohseni' via deal.II User Group]

Thanks for your suggestions.


Generally, in FEM you number elements like in the below picture: 




I understand of course, that the cells (elements) in deal.II are parallely 
distributed, but it should be possible to identify the geometrical position 
and ID of the elements and their corresponding nodes. 


For instance the following works:


typename parallel::distributed::Triangulation::active_cell_iterator 
cell = triangulation.begin_active(), endc = triangulation.end(); 
for (; cell != endc; ++cell)
{
 std::cout << "CELL NUMBER: " << cell->id() << std::endl;
}



Of course here the locally owned term is missing, but is it correct to 
access the cells like this?


I am a bit confused since it is a rather easy task which should not be such 
a big "deal" in deal.II.


Best 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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Access specific element within a distributed triangulation

['Seyed Ali Mohseni' via deal.II User Group]

Hi,

Is it possible to access a specific element within a parallel distributed 
triangulation by means of the following procedure for instance:

typename DoFHandler::active_cell_iterator cell = 
dof_handler.begin_active(), endc = dof_handler.end();

for (; cell != endc; ++cell)
{
 if (cell->is_locally_owned())
 {
 std::cout << "CELL NUMBER 4 VERTEX 1: " << cell(4)->vertex(0) << std::endl;
 }
}

The above won't compile.

How am I able to access for example element number 4 in order to obtain the 
coordinates of point 1 within the element.

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Synchronization issue for parallel computation of material forces

['Seyed Ali Mohseni' via deal.II User Group]

Maybe this helps: 

How do you compute "reaction forces" within deal.II? 

Did you ever compute reaction forces and write them as output within a 
parallel framework?

Best,
Seyed Ali

-- 
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.

[deal.II] Re: Synchronization issue for parallel computation of material forces

['Seyed Ali Mohseni' via deal.II User Group]

Thanks, Jean-Paul. I was thinking maybe something in the same way like 
step-40.

I tried your approach, but can I just set it equal to each other?

I did the following:

*Main class*
LA::MPI::Vector configurational_forces;
LA::MPI::Vector local_configurational_forces;

*setup_system() class*
configurational_forces.reinit(locally_owned_dofs, locally_relevant_dofs, 
mpi_com);

*compute_config_forces() class*
Vector cell_cf(dofs_per_cell); // Cell Configurational forces

typename DoFHandler::active_cell_iterator cell = 
dof_handler.begin_active(), endc = dof_handler.end();
for (; cell != endc; ++cell)
{
if ( cell->is_locally_owned() )
{
fe_values.reinit(cell);

cell_cf = 0;

...

// Configurational force computation ...

...

cell->get_dof_indices(local_dof_indices);

constraints.distribute_local_to_global(cell_cf, 
local_dof_indices, *local_configurational_forces*);
}
}

*local_configurational_forces*.compress(VectorOperation::add);

*output_results() class*
// [ CONFIGURATIONAL FORCES 
]==

configurational_forces = local_configurational_forces;  *(here I receive an 
error!)*

std::vector configurational_forces_magnitude(dim, 
"config_forces");

std::vector 
configurational_forces_interpretation(dim, 
DataComponentInterpretation::component_is_part_of_vector);

data_out.add_data_vector(configurational_forces, 
configurational_forces_magnitude, DataOut::type_dof_data, 
configurational_forces_interpretation);


I receive the following error:

make all 
Scanning dependencies of target solid_mechanics
[ 50%] Building CXX object 
CMakeFiles/solid_mechanics.dir/solid_mechanics.cc.o
/home/seyedali/fe_models/deal.II/solid_mechanics/solid_mechanics.cc(950): 
error: no operator "=" matches these operands
operand types are: const 
dealii::LinearAlgebraPETSc::MPI::Vector = const 
dealii::LinearAlgebraPETSc::MPI::Vector
  configurational_forces = local_configurational_forces;
^
  detected during:
instantiation of "void 
SolidMechanics::do_initial_timestep() [with dim=2]" at line 1081
instantiation of "void SolidMechanics::run() [with dim=2]" 
at line 1113

compilation aborted for 
/home/seyedali/fe_models/deal.II/solid_mechanics/solid_mechanics.cc (code 2)
make[2]: *** [CMakeFiles/solid_mechanics.dir/solid_mechanics.cc.o] Fehler 2
CMakeFiles/solid_mechanics.dir/build.make:62: die Regel für Ziel 
„CMakeFiles/solid_mechanics.dir/solid_mechanics.cc.o“ scheiterte
CMakeFiles/Makefile2:195: die Regel für Ziel 
„CMakeFiles/solid_mechanics.dir/all“ scheiterte
make[1]: *** [CMakeFiles/solid_mechanics.dir/all] Fehler 2
Makefile:83: die Regel für Ziel „all“ scheiterte
make: *** [all] Fehler 2

I assume the size doesn't match. Can I really just write it like you 
suggested?

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Synchronization issue for parallel computation of material forces

['Seyed Ali Mohseni' via deal.II User Group]

And this is the procedure how I output them:

DataOut data_out;
data_out.attach_dof_handler(dof_handler);

std::vector configurational_forces_magnitude(dim, 
"config_forces");

std::vector 
configurational_forces_interpretation(dim, 
DataComponentInterpretation::component_is_part_of_vector);

data_out.add_data_vector(configurational_forces, 
configurational_forces_magnitude, DataOut::type_dof_data, 
configurational_forces_interpretation);

Very weird. I can compute everything correctly for one core, but why is it 
not multicore compatible?


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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Synchronization issue for parallel computation of material forces

['Seyed Ali Mohseni' via deal.II User Group]

Hi,

With regard to the approach which was recently discussed in: 
https://groups.google.com/d/msg/dealii/25C57MYvfpg/JFB8aXKoAgAJ

There is a little problem when computing the material forces in a parallel 
setting.


  
  



I assume it has something to do with the synchronization. 

The computation of material forces is done similar to the system_rhs vector.



Vector cell_cf(dofs_per_cell); // Cell Configurational forces

typename DoFHandler::active_cell_iterator cell = 
dof_handler.begin_active(), endc = dof_handler.end();
for (; cell != endc; ++cell)
{
if ( cell->is_locally_owned() )
{
fe_values.reinit(cell);

cell_cf = 0;

...

// Configurational force computation ...

...

cell->get_dof_indices(local_dof_indices);

constraints.distribute_local_to_global(cell_cf, 
local_dof_indices, configurational_forces);
}
}

configurational_forces.compress(VectorOperation::add);


For a single core computation it is validated and correct.


What could be possibly wrong when going to parallel?



Best regards,

Seyed Ali 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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Shape function derivatives

['Seyed Ali Mohseni' via deal.II User Group]

@Prof. Bangerth: Thank you very much :) Finally, I understand how the 
stiffness matrix in step-18 is constructed. Your shape function derivatives 
are basically already the B-operator. 

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Shape function derivatives

['Seyed Ali Mohseni' via deal.II User Group]

Hi,

I found some issue when constructing the FE B-operator for multiple 
elements. In a single element it works perfectly.
My approach is the following:

template
FullMatrix get_b_operator (const FEValues _values, const 
unsigned int dofs_per_cell, const unsigned int q)
{
FullMatrix tmp(dim, GeometryInfo::vertices_per_cell);

std::vector > invJ = 
fe_values.get_inverse_jacobians();

for (unsigned int i = 0; i < dofs_per_cell; i = i + dim)
{
const unsigned int index = i / dim;

// COMPUTE: B-operator (Remark: This version has to be extended for 3D!)
tmp[0][index] = invJ[q][0][0] * fe_values.shape_grad(i, q)[0] + 
invJ[q][0][1] * fe_values.shape_grad(i, q)[1];
tmp[1][index] = invJ[q][1][0] * fe_values.shape_grad(i, q)[0] + 
invJ[q][1][1] * fe_values.shape_grad(i, q)[1];
}

return tmp;
}

Hence, I am computing B = inv(J) * transpose(dN). For a single element I 
receive the following correct results for dN:

*C++ CODE*


SHAPE FUNCTION DERIVATIVES (GAUSS POINT 1)

-0.788675  -0.788675  
0.788675  -0.211325  
0.211325  0.211325  
-0.211325  0.788675  

SHAPE FUNCTION DERIVATIVES (GAUSS POINT 2)

-0.788675  -0.211325  
0.788675  -0.788675  
0.211325  0.788675  
-0.211325  0.211325  

SHAPE FUNCTION DERIVATIVES (GAUSS POINT 3)

-0.211325  -0.788675  
0.211325  -0.211325  
0.788675  0.211325  
-0.788675  0.788675  

SHAPE FUNCTION DERIVATIVES (GAUSS POINT 4)
-0.211325  -0.211325  
0.211325  -0.788675  
0.788675  0.788675  
-0.788675  0.211325  



*DEAL.II*


SHAPE FUNCTION DERIVATIVES

  -0.788675 -0.788675   -0.788675 -0.788675   0.788675 -0.211325   0.788675 
-0.211325   -0.211325 0.788675   -0.211325 0.788675   0.211325 0.211325   
0.211325 0.211325 
  -0.788675 -0.211325   -0.788675 -0.211325   0.788675 -0.788675   0.788675 
-0.788675   -0.211325 0.211325   -0.211325 0.211325   0.211325 0.788675   
0.211325 0.788675 
  -0.211325 -0.788675   -0.211325 -0.788675   0.211325 -0.211325   0.211325 
-0.211325   -0.788675 0.788675   -0.788675 0.788675   0.788675 0.211325   
0.788675 0.211325 
  -0.211325 -0.211325   -0.211325 -0.211325   0.211325 -0.788675   0.211325 
-0.788675   -0.788675 0.211325   -0.788675 0.211325   0.788675 0.788675   
0.788675 0.788675 


Due to the ordering of the element connectivity (numbering of nodes are 
different in deal.II) we have to swap column 3 and 4 for comparison.

Now for 4 elements there is something strange when I compute my shape 
function derivatives within deal.II:

*DEAL.II*


SHAPE FUNCTION DERIVATIVES

  -1.577350 -1.577350   -1.577350 -1.577350   1.577350 -0.422650   1.577350 
-0.422650   -0.422650 1.577350   -0.422650 1.577350   0.422650 0.422650   
0.422650 0.422650 
  -1.577350 -0.422650   -1.577350 -0.422650   1.577350 -1.577350   1.577350 
-1.577350   -0.422650 0.422650   -0.422650 0.422650   0.422650 1.577350   
0.422650 1.577350 
  -0.422650 -1.577350   -0.422650 -1.577350   0.422650 -0.422650   0.422650 
-0.422650   -1.577350 1.577350   -1.577350 1.577350   1.577350 0.422650   
1.577350 0.422650 
  -0.422650 -0.422650   -0.422650 -0.422650   0.422650 -1.577350   0.422650 
-1.577350   -1.577350 0.422650   -1.577350 0.422650   1.577350 1.577350   
1.577350 1.577350 



 The results for one element change now by a factor of 2 and with 
increasing element size it always differs by another factor of 2 again.

I was wondering, if my approach to compute my B-operator was wrong. Doesn't 
fe_values.shape_grad() give me the shape function derivatives directly? Or 
do I have to use shape_grad_component? 

Additionally, I checked the symmetric_gradient from step-18 and it contains 
actually the correct values I need, but somehow in a strange order:


SYMMETRIC GRADIENT (GAUSS POINT 1)

0.00 0.211325 
0.211325 0.422650

SYMMETRIC GRADIENT (GAUSS POINT 2)

0.00 0.211325 
0.211325 1.577350

SYMMETRIC GRADIENT (GAUSS POINT 3)

0.00 0.788675 
0.788675 0.422650 

SYMMETRIC GRADIENT (GAUSS POINT 4)

0.00 0.788675 
0.788675 1.577350




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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Assembly of material forces

['Seyed Ali Mohseni' via deal.II User Group]

@Prof. Bangerth: The upper one where pseudocolor is written, is from VisIt 
and deal.II. Below as a reference I showed my results from the C++ code. 
Both simulations are done using Q4 quadrilateral elements. There are no 
triangular elements here and the meshes are identical. They are 
intentionally created for comparison.

Here I show a more refined version, maybe you see the solution more clearly:






Additionally, I colored the configurational force arrows from the deal.II 
solution to black, so it becomes more visible.

I think I found the last small detail which seemed to be wrong, but am 
still not 100 % sure. There is something I don't understand completely. 
Surprisingly, the JxW value in deal.II stays the same value like in my C++ 
code, although the jacobian matrix differs by a factor of 2. Exactly this 
factor of 2 was the difference in my previous results. Now as you can see I 
get identical results in both approaches. 

I also fixed some other problem, namely the storage of my data for each 
element separately. Earlier I stored everything in Gauss points only, but 
the global vector for the PointHistory was overwritten since I had now 
several elements. For a single element it worked like a charm which is why 
I didn't recognize it directly. Now I am curious, if this is the best way 
to store data in a vector for large number of elements. I think you have 
probably some deal.II way of storing data on cells or at least efficiently 
in a vector using pointer arithmetics. Is the CellStorageData Jean-Paul 
mentioned what I seek? I really don't want to lose the scalability of my 
approach since it should be applicable to large problems also.


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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Assembly of material forces

['Seyed Ali Mohseni' via deal.II User Group]

@Prof. Bangerth: The upper one where pseudocolor is written, is from VisIt 
and deal.II. Below as a reference I showed my results from the C++ code. 
Both simulations are done using Q4 quadrilateral elements. There are no 
triangular elements here and the meshes are identical. They are 
intentionally created for comparison.

Here I show a more refined version, maybe you see the solution more clearly:






Additionally, I colored the configurational force arrows from the deal.II 
solution to black, so it becomes more visible.

I think I found the last small detail which seemed to be wrong, but am 
still not 100 % sure. There is something I don't understand completely. 
Surprisingly, the JxW value in deal.II stays the same value like in my C++ 
code, although they jacobian matrix differs by a factor of 2. Exactly this 
factor of 2 was the difference in my prior results. Now as you can see I 
get identical results in both approaches. 

I also fixed some other problem, namely the storage of my necessary data 
for each element separately. Earlier I stored everything in Gauss points 
only, but the global vector was overwritten since I had now several 
elements. For a single element it worked like a charm which is why I didn't 
recognize it directly. Now I am curious, if this is the best way to store 
data in a vector for large number of elements. I think you have probably 
some deal.II way of storing data on cells or at least efficiently in a 
vector using pointer arithmetics. Is the CellStorageData Jean-Paul 
mentioned what I seek? I really don't want to lose the scalability of my 
approach since it should be applicable to large problems also.


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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Assembly of material forces

['Seyed Ali Mohseni' via deal.II User Group]

First of all, thank you for all your suggestions. They helped me a lot, 
especially from Prof. Bangerth and Jean-Paul. Many thanks ;)

I solved the assembly by listening to Prof. Bangerth's advice that I should 
try to implement things in a deal.II way. And after Jean-Paul told me to 
take a look at step-42 I realized how things are working in deal.II.
Unfortunately, I have still some issues with the values. They are not the 
same. I believe this is dependent on the B-operator, namely the jacobians. 
They differ by factor 2 probably since our unit cell in deal.II is within 
[0,1]^dim and the C++ code uses [-1,1]^dim.
But shouldn't the assembled configurational forces still be the same 
independent from the mapping?

Here is what I got so far :)





Kind regards,
Seyed Ali

-- 
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.

[deal.II] Re: Assembly of material forces

['Seyed Ali Mohseni' via deal.II User Group]

First of all, thank you for all your suggestions. They helped me a lot, 
especially from Prof. Bangerth and Jean-Paul. Many thanks ;)

I somehow figured the assembly out by listening to Prof. Bangerth's advice 
that I should try to solve it the deal.II way. And after Jean-Paul told me 
to take a look at step-42 I realized how things are working in deal.II.
Unfortunately, I have still some issues with the values. They are not the 
same. I believe this is dependent on the B-operator, namely the jacobians. 
They differ by factor 2 probably since our unit cell in deal.II is within 
[0,1]^dim and the C++ code uses [-1,1]^dim.
But shouldn't the assembled configurational forces still be the same 
independent from the mapping?

Here is what I got so far :)





Kind regards,
Seyed Ali

-- 
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.

[deal.II] Re: Assembly of material forces

['Seyed Ali Mohseni' via deal.II User Group]

@Daniel:

So you have a FESystem<2>(FE_Q<2>(1),2) element, i.e. linear continuous 
> elements with two components, and the vector represents corresponding local 
> DoF values?
>

 Exactly. I defined first 

FESystem fe;

and then 

template
SolidMechanics::SolidMechanics(): ..., fe(FE_Q(1), dim), ...

as constructor initializer. 
Here dim=2 would result in what you just mentioned.

In particular, assigning the value to the DoF with the corresponding 
> support point and component is all you have to do.
>

Is there a command that allows me to assign values to the DoF. Can you 
please explain this a bit more.

The Gauss quadrature points never include the vertices [1] so saying that 
> the values represent forces for each vertex and at the same time that the 
> values are computed within Gauss points sounds contradictory.
>

Yes, you are right. I described it a bit mixed up. The vector is just 
shaped and designed for our global DoFs structure. You initialize it by:

 
Vector cell_cf(dofs_per_cell);

Of course integration points and vertices are different. Thanks for the 
reference, it may come useful.

Note that this approach is only suitable for Q1 elements. For higher 
polynomial degree, you have to calculate the values for all the support 
points of the finite element.

Would you suggest another approach? Since it may be useful to implement 
this in a general way, useable for higher degrees of freedom. 


@Jean-Paul:

Providing the theory is for our benefit (and therefore yours). No one here 
> is necessarily an expert in what you're trying to accomplish. So, for all 
> of your explanation, an equation or two might go a long way to help us 
> understand exactly what you're trying to achieve and, therefore, help us 
> suggest how you would accomplish it. If you can't express your problem 
> clearly then please don't be surprised if you don't receive any meaningful 
> assistance. Help us help you.
>

Of course, but I just wanted to save you from additional headaches :) 
The computation of configurational forces is accomplished by the following 
formula:




Deal.II does not attempt to be the beginning and end of provided 
> functionality. Yes, there are times when deal.II does not provide some 
> functionality that you require. There are a zoo of other libraries that you 
> could also use to extend your code and fill in the gaps that deal.II does 
> not provide. Of course, we always welcome extensions to existing features 
> (e.g the relatively new CellDataStorage class 
> 
>  and TransferableQuadraturePointData class 
> )
>  
> or the addition of new ones.
>

How can I use these new features? Is there a step example or tutorial? Is 
this TransferableQuadraturePointData the same like in step-18?

 Is this a deal.II error or a user error? Have you repopulate the 
> local_dof_indices vector for each new cell that you're doing this 
> "assembly" on? It looks like all contributions are going to a fixed set of 
> global DoFs.


I think more of a user error. Deal.II is nice and very elegantly designed. 
I just don't think it will automatically assemble it correctly in this 
case. First I have to understand it, before telling the program. ;)

I just did this:

std::vector local_dof_indices(dofs_per_cell);

typename DoFHandler::active_cell_iterator cell = dof_handler.
begin_active(), endc = dof_handler.end();
for (; cell != endc; ++cell)
{
if ( cell->is_locally_owned() )
{
   fe_values.reinit(cell);

 ... configurational forces computation ...

}

configurational_forces = cell_cf;

cell->get_dof_indices(local_dof_indices);

constraints.distribute_local_to_global(cell_cf, local_dof_indices, 
configurational_forces);
}


I am really not sure, if the hanging node constraints have anything to do 
with my desired assembly at the moment. I assume they play a role when I 
refine the mesh.
So what exactly do you mean by repopulating each cell?

Many thanks, this group here is full of nice people. Hopefully, I can help 
you also as a reward in future. :)

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Assembly of material forces

['Seyed Ali Mohseni' via deal.II User Group]

@Daniel:

So you have a FESystem<2>(FE_Q<2>(1),2) element, i.e. linear continuous 
> elements with two components, and the vector represents corresponding local 
> DoF values?
>

 Exactly. I defined first 

FESystem fe;

and then 

template
SolidMechanics::SolidMechanics(): ..., fe(FE_Q(1), dim), ...

as constructor initializer. 
Here dim=2 would result in what you just mentioned.

In particular, assigning the value to the DoF with the corresponding 
> support point and component is all you have to do.
>

Is there a command that allows me to assign values to the DoF. Can you 
please explain this a bit more.

The Gauss quadrature points never include the vertices [1] so saying that 
> the values represent forces for each vertex and at the same time that the 
> values are computed within Gauss points sounds contradictory.
>

Yes, you are right. I described it a bit mixed up. The vector is just 
shaped and designed for our global DoFs structure. You initialize it by:

 
Vector cell_cf(dofs_per_cell);

Of course integration points and vertices are different. Thanks for the 
reference, it may come useful.

Note that this approach is only suitable for Q1 elements. For higher 
polynomial degree, you have to calculate the values for all the support 
points of the finite element.

Would you suggest another approach? Since it may be useful to implement 
this in a general way, useable for higher degrees of freedom. 


@Jean-Paul:

Providing the theory is for our benefit (and therefore yours). No one here 
> is necessarily an expert in what you're trying to accomplish. So, for all 
> of your explanation, an equation or two might go a long way to help us 
> understand exactly what you're trying to achieve and, therefore, help us 
> suggest how you would accomplish it. If you can't express your problem 
> clearly then please don't be surprised if you don't receive any meaningful 
> assistance. Help us help you.
>

Of course, but I just wanted to save you from additional headaches :) 
The computation of configurational forces is accomplished by the following 
formula:




Deal.II does not attempt to be the beginning and end of provided 
> functionality. Yes, there are times when deal.II does not provide some 
> functionality that you require. There are a zoo of other libraries that you 
> could also use to extend your code and fill in the gaps that deal.II does 
> not provide. Of course, we always welcome extensions to existing features 
> (e.g the relatively new CellDataStorage class 
> 
>  and TransferableQuadraturePointData class 
> )
>  
> or the addition of new ones.
>

How can I use these new features? Is there a step example or tutorial? Is 
this TransferableQuadraturePointData the same like in step-18?

 Is this a deal.II error or a user error? Have you repopulate the 
> local_dof_indices vector for each new cell that you're doing this 
> "assembly" on? It looks like all contributions are going to a fixed set of 
> global DoFs.


I think more of a user error. Deal.II is nice and very elegantly designed. 
I just don't think it will automatically assemble it correctly in this 
case. First I have to understand it, before telling the program. ;)

I just did this:

std::vector local_dof_indices(dofs_per_cell);

typename DoFHandler::active_cell_iterator cell = 
dof_handler.begin_active(), endc = dof_handler.end();
for (; cell != endc; ++cell)
{
 if ( cell->is_locally_owned() )
 {

... configurational forces computation ...



 

-- 
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.

[deal.II] Re: Assembly of material forces

['Seyed Ali Mohseni' via deal.II User Group]

My problem is not the theory. I am just not a deal.II expert and therefore 
it is still quite difficult for me to implement my approach. Going in to 
theoretical detail would take maybe much longer. Maybe it is better 
instead, if I reformulate my question like this:

There is a vector called cell_cf for instance. It contains all integrated 
material forces (summed up material forces from integration points) for 
each cell. Each cell has 4 vertices and 8 DoFs in 2D. Hence, it has to be a 
8x1 vector.

For instance:

CONFIGURATIONAL FORCES WITHIN ELEMENT
   0.010126  
   0.007502  
   0.011509  
   -0.012249  
   -0.008940  
   -0.009762  
   -0.012695  
   0.014508 

Now imagine you have this vector which represents forces for each vertex, 
e.g. vertex 1 has 0.010126 in X-direction and 0.007502 in Y-direction, 
vertex 2 has  0.011509 in X-direction and -0.012249  in Y-direction and so 
on...

How is one able now to take these vectors for each cell and assemble them 
correspondent to the DoF and summarize the values for coinciding nodes. 
Remember the values are computed within Gauss points and are just 
summarized over integration points to one single vector for each cell. It 
means there should be nothing stored on the nodes yet, but we have to take 
the information to the global structure and node/vertices level.

I already tried to take the approach we implemented within our C++ code, 
but it is quite difficult to store data in deal.II since there is no 
processinfo or such a thing which tracks all data created and stored over 
the computation. I cannot save my configurational force vector and access 
it anywhere except I store it by using the setup_quadrature_point_history() 
 and update_quadrature_point_history() functions from step-18. But this 
means I have to store my elastic energy then update and compute my material 
forces which are based on the energy and update again. This is not very 
nice to call a function several times for every little thing you want to 
save.

Isn't there a better way to store data efficiently during runtime?

Sorry for causing headaches ^^'

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Assembly of material forces

['Seyed Ali Mohseni' via deal.II User Group]

Dear Prof. Bangerth,

It fails, because I tried to use the same procedure for the standard 
assembly such as in step-18 or step-40 and I received some awkward contour 
plots which won't match with what I showed from the C++ code. Hence, fails 
is refered to my own failure in understanding the core structure of deal.II 
;)

Can I assemble a vector such as these configurational force vector (8x1) 
after I solved the system? Because my compute_config_forces() function is 
called after I solved the system. Hence, I am wondering about a procedure 
such as:

configurational_forces = cell_cf;

cell->get_dof_indices(local_dof_indices);

constraints.distribute_local_to_global(cell_cf, local_dof_indices, 
configurational_forces);


where configurational_forces represents a Vector variable I defined 
globally within my class.

I also found some MeshWorker approach such as shown in step-12. Is it 
possible to use such a method to reach my aim?

Conclusion: It would be great, if you can at least guide me to the correct 
direction of how to assemble a vector after I solved my system. Assembly of 
the cell vector should also take all existing boundary conditions, 
constraints etc. into account.

Thank you :)

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Assembly of material forces

['Seyed Ali Mohseni' via deal.II User Group]

Dear deal.II experts,

I have a question regarding the assembly of material forces also known as 
configurational forces. 

First of all, I have constructed my own solid mechanics code within deal.II 
which I already validated with our own C++ code. Hence, the deal.II code 
works until here.








Attached you find my results on material forces for one single element. It 
also works perfectly. 




Now I try to compute the material forces for the whole structure. Hence, I 
tried using the same procedure of assembly for our LHS or RHS matrices in 
deal.II, but it fails. Do I have to take a special approach in order to 
make it work? I assume my local material forces have to be assembled to the 
global nodes, but this assembly is it just a simple mapping to the nodes? 
Since there are quadtree constraints, it will not be that simple, I think.





The successful results are computed within our simple c++ code and I just 
show them for comparison reasons. 

Extra question: Would my work be of interest for your code gallery maybe in 
the near future?  


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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Postprocessing using the same variable twice

['Seyed Ali Mohseni' via deal.II User Group]

Thank you, Daniel.

I solved my problem already in the following way before reading your post:

DataOut data_out;
data_out.attach_dof_handler(dof_handler);

// [ DISPLACEMENT COMPONENTS 
]=

std::vector displacement_components;
switch (dim)
{
case 1:
displacement_components.push_back("displacement_x");
break;
case 2:
displacement_components.push_back("displacement_x");
displacement_components.push_back("displacement_y");
break;
case 3:
displacement_components.push_back("displacement_x");
displacement_components.push_back("displacement_y");
displacement_components.push_back("displacement_z");
break;
default:
Assert(false, ExcInternalError())
;
}

std::vector 
displacement_component_interpretation(dim, 
DataComponentInterpretation::component_is_scalar);

// Collect data output
data_out.add_data_vector(incremental_displacement, displacement_components, 
DataOut::type_dof_data, displacement_component_interpretation);

// [ DISPLACEMENTS MAGNITUDE 
]=

std::vector displacement_magnitude(dim, "displacement");

std::vector 
displacement_magnitude_interpretation(dim, 
DataComponentInterpretation::component_is_part_of_vector);

data_out.add_data_vector(incremental_displacement, displacement_magnitude, 
DataOut::type_dof_data, displacement_magnitude_interpretation);


I thought there would be a more elegant solution than this, which is why I 
asked here. But it seems your suggestion is similar to my approach. 
Regarding the DataPostprocessor class, I am just using it for quantities at 
integration points, e.g. strains etc. 

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Postprocessing using the same variable twice

['Seyed Ali Mohseni' via deal.II User Group]

Dear all,

Is it possible to output the magnitude of displacements and its components 
all together?

std::vector 
data_component_interpretation(dim, DataComponentInterpretation::
component_is_part_of_vector);

The above allows me to output the displacements as a vector field and shows 
one variable containing the magnitude of displacements.

std::vector 
data_component_interpretation(dim, DataComponentInterpretation::
component_is_scalar);

Hence, this one enables the output of each component separately, e.g. 
displacement_x (delta_x) according to step-18.

What I want is to have all in one, such as:

displacement_x
displacement_y
displacement_z
displacement_magnitude

Is that possible?

I am using VisIt 2.11.0.

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Elastic energy computation

['Seyed Ali Mohseni' via deal.II User Group]

Thanks, I understand now. 

I just deleted the post, because I misunderstood what you meant earlier. 
I thought you are refering to the left eps_phi_i (first) and right 
eps_phi_j (second).

But you were talking about: 

template  inline SymmetricTensor<2,dim> get_strain (const 
std::vector > )

This function differs of course since it computes the strain based on the 
displacement gradients. This version I already computed, as I mentioned 
before.
I thought maybe you have some procedure to store the last converged 
displacements so we can use them for the computation of strains within the 
assembly routine.
This seems not easily possible due to the structure of deal.II. From what I 
understood we have to solve first before obtaining the "physical" strains.

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Elastic energy computation

['Seyed Ali Mohseni' via deal.II User Group]

I checked the second one. It doesn't coincide with the strain tensor I 
obtain from the energy computation after I solve the system, e.g.

 - in GAUSS point 1 - 

ASSEMBLY STRAIN TENSOR
0 0.105662 
0.105662 0.211325 

ENERGY STRAIN TENSOR
 0.001081   -0.000738  
 -0.000738   -0.008525 

... etc.

How do you assume that it is the same strain tensor?

The PointHistory idea is very good. Thanks for remembering me ;)

Best regards,
Seyed Ali
 

-- 
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.

[deal.II] Elastic energy computation

['Seyed Ali Mohseni' via deal.II User Group]

Hi,

Is it possible to compute the elastic energy for a simple 2D cube before we 
solve the system?
Because I checked the displacement gradients and they are 0 due to having 
no solution probably.

Hence, I wonder how would someone be able to compute the elastic energy 
with the following formula:

psi = 0.5 * lambda * tr_eps * tr_eps + mu * tr_eps_2;

with 

tr_eps = trace(eps);
tr_eps_2 = trace(eps * eps);

Here eps is the strain tensor and lambda as well as mu are the lamé 
parameters, respectively.

In order to compute the strain tensor eps I need the displacement gradients:

 eps = 0.5 * (grad_u + transpose(grad_u));

And I from what I checked, the get_strain() function from step-18 doesn't 
give me the same strain tensor. It seems more like a B-operator to me.

If I compute the elastic energy after I solved the system, how am I able to 
store the energy values for each cell, so I can have access to it all the 
time. 

>From Timo Heister and Thomas Wick's phase-field crack examples I learned 
that they compute energy by means of the compute_energy() function, but it 
seems to me like they just output the value and it is not used further 
during computation.

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Strain tensor and displacement gradients

['Seyed Ali Mohseni' via deal.II User Group]

I checked the Gauss points within my normal computation using FEValues and 
compared them with the evaluation points from the DataPostprocessor:

GAUSS POINT COORDINATES

  0.211325   0.211325 
  0.788675   0.211325 
  0.211325   0.788675 
  0.788675   0.788675

EVALUATION POINTS

  0.00   0.00 
  1.00   0.00 
  0.00   1.00 
  1.00   1.00

This means everything we do in DataPostprocessor is evaluated directly in 
the mapped Gauss points which is actually the reason why we use this 
postprocessor approach. We want to map the solution of quantities stored on 
the Gauss points to the physical nodes, so we can visualize them in a 
proper way.

I hope this is correct, because then I assume that both approaches are 
correct and just differ by a simple mapping procedure.

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Strain tensor and displacement gradients

['Seyed Ali Mohseni' via deal.II User Group]

>From what I learned, all information we obtain from FEValues is 
automatically mapped from the unit cell to the physical space/coordinates. 
By patch vertices you mean the nodes of each element, hence 
DataPostprocessor gives me the gradients with respect to the mapped gauss 
points to the nodes of each element? And standard FEValues gives me the 
gradients with respect to the mapped physical coordinates. Isn't this 
somehow the same or where is the key difference here? Sorry for being 
stupid right now ^^'

Best regards,
Seyed Ali

-- 
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.

[deal.II] Re: Strain tensor and displacement gradients

['Seyed Ali Mohseni' via deal.II User Group]

>From what I learned, all informtion we obtain from FEValues is 
aoutomatically mapped from the unit cell to the physical space/coordinates. 
By patch vertices you mean the nodes of each element, hence 
DataPostprocessor gives me the gradients with respect to the mapped gauss 
points to the nodes of each element? And standard FEValues gives me the 
gradients with respect to the mapped physical coordinates. Isn't this 
somehow the same or where is the key difference here? Sorry for being 
stupid right now ^^'

Best regards,
Seyed Ali

-- 
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.

[deal.II] Re: Examples makefile

['Seyed Ali Mohseni' via deal.II User Group]

Hi Victor,

I would check first of all, if you compiled Trilinos with the correct 
settings, e.g.

mkdir build
cd build

cmake \
-DCMAKE_C_COMPILER=/path_to_your_installation/mpich-3.2/bin/mpicc   
\
-DCMAKE_CXX_COMPILER=/path_to_your_installation/mpich-3.2/bin/mpicxx 
\
-DCMAKE_Fortran_COMPILER=/path_to_your_installation/mpich-3.2/bin/mpifort   
\
-DTrilinos_ENABLE_Amesos=ON   \
-DTrilinos_ENABLE_Epetra=ON   \
-DTrilinos_ENABLE_Ifpack=ON   \
-DTrilinos_ENABLE_AztecOO=ON \
-DTrilinos_ENABLE_Sacado=ON   \
-DTrilinos_ENABLE_Teuchos=ON \
-DTrilinos_ENABLE_MueLu=ON   \
-DTrilinos_ENABLE_ML=ON   \
-DTrilinos_VERBOSE_CONFIGURE=OFF \
-DTrilinos_ENABLE_EXPLICIT_INSTANTIATION=ON   \
-DTPL_ENABLE_MPI=ON   \
-DBUILD_SHARED_LIBS=ON   \
-DCMAKE_VERBOSE_MAKEFILE=OFF \
-DCMAKE_BUILD_TYPE=RELEASE   \
-D 
CMAKE_INSTALL_PREFIX:PATH=/path_to_the_folder_where_trilinos_should_be_installed/trilinos-12.8.1
 
\
../

make -j8 install

Then I would see, if python is configured correctly. Maybe your python 
libraries are buggy, so better check them with other programs that use 
python.

Another option, would be to use candi: https://github.com/koecher/candi

Hope these advices may be of little help to you.

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Strain tensor and displacement gradients

['Seyed Ali Mohseni' via deal.II User Group]

Hi,

Using the postprocessor approach according to step-33 I have computed the 
strain tensor successfully by means of the 
Postprocessor::compute_derived_quantities_vector function and its 
displacement gradient duh.

Unfortunately, within an arbitrary function called right after I solve the 
system, exactly before postprocessing, I try the same, but obtain different 
results:

const QGauss quadrature_formula(2);

const unsigned int n_q_points = quadrature_formula.size();

 FEValues fe_values(fe, quadrature_formula, update_values |
 update_gradients | update_quadrature_points | update_JxW_values);

 LA::MPI::Vector rel_solution;
rel_solution = locally_relevant_solution;

 std::vector solution_values(n_q_points, Vector(
dim + 1));

 std::vector > > solution_grads(n_q_points, std::
vector >(dim + 1));

 typename DoFHandler::active_cell_iterator cell = dof_handler.
begin_active(), endc = dof_handler.end();
for (; cell != endc; ++cell)
{
  if ( cell->is_locally_owned() )
  {
   fe_values.reinit(cell);

fe_values.get_function_values(rel_solution, solution_values);
   fe_values.get_function_gradients(rel_solution, solution_grads);

Tensor<2, dim> grad_u, eps;

for (unsigned int q = 0; q < n_q_points; ++q)
   {
 std::cout << "\n - in GAUSS point " << q + 1 << " - \n" <<
 std::endl;

  for (unsigned int i = 0; i < dim; ++i)
 {
   grad_u[i] = solution_grads[q][i];

eps = 0.5 * (grad_u + transpose(grad_u));
  }

  print_tensor(grad_u, "DISPLACEMENT GRADIENT");
  print_tensor(eps, "STRAIN TENSOR");
  }

  }
}


The above code gives me the following output:


- in GAUSS point 1 - 

DISPLACEMENT GRADIENT
 0.001081   -0.001477  
 -0.00   -0.008525  

STRAIN TENSOR
 0.001081   -0.000738  
 -0.000738   -0.008525  


 - in GAUSS point 2 - 

DISPLACEMENT GRADIENT
 0.001081   0.001477  
 -0.00   -0.008525  

STRAIN TENSOR
 0.001081   0.000738  
 0.000738   -0.008525  


 - in GAUSS point 3 - 

DISPLACEMENT GRADIENT
 0.004034   -0.001477  
 -0.00   -0.008525  

STRAIN TENSOR
 0.004034   -0.000738  
 -0.000738   -0.008525  


 - in GAUSS point 4 - 

DISPLACEMENT GRADIENT
 0.004034   0.001477  
 -0.00   -0.008525  

STRAIN TENSOR
 0.004034   0.000738  
 0.000738   -0.008525  


And the below are results from the postprocessor computation for a single 
element (already checked them, they should be correct):


- in GAUSS point 1 - 

DISPLACEMENT GRADIENT
 0.00   -0.002557  
 0.00   -0.008525  

STRAIN TENSOR
 0.00   -0.001279  
 -0.001279   -0.008525  


 - in GAUSS point 2 - 

DISPLACEMENT GRADIENT
 0.00   0.002557  
 0.00   -0.008525  

STRAIN TENSOR
 0.00   0.001279  
 0.001279   -0.008525  


 - in GAUSS point 3 - 

DISPLACEMENT GRADIENT
 0.005115   -0.002557  
 -0.00   -0.008525  

STRAIN TENSOR
 0.005115   -0.001279  
 -0.001279   -0.008525  


 - in GAUSS point 4 - 

DISPLACEMENT GRADIENT
 0.005115   0.002557  
 -0.00   -0.008525  

STRAIN TENSOR
 0.005115   0.001279  
 0.001279   -0.008525  

Why is the fe_values.get_function_gradients() function giving me wrong 
gradients?

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Strain tensor and displacement gradients

['Seyed Ali Mohseni' via deal.II User Group]

Hi,

Using the postprocessor approach according to step-33 I have computed the 
strain tensor successfully by means of the 
Postprocessor::compute_derived_quantities_vector function and its 
displacement gradient duh.

Unfortunately, within an arbitrary function called right after I solve the 
system, exactly before postprocessing, I try the same, but obtain different 
results:

const QGauss quadrature_formula(2);

const unsigned int n_q_points = quadrature_formula.size();

 FEValues fe_values(fe, quadrature_formula, update_values | 
update_gradients | update_quadrature_points | update_JxW_values);

 LA::MPI::Vector rel_solution;
rel_solution = locally_relevant_solution;

 std::vector solution_values(n_q_points, Vector(dim 
+ 1));

 std::vector > > solution_grads(n_q_points, std::
vector >(dim + 1));

 typename DoFHandler::active_cell_iterator cell = dof_handler.
begin_active(), endc = dof_handler.end();
for (; cell != endc; ++cell)
{
  if ( cell->is_locally_owned() )
  {
   fe_values.reinit(cell);

fe_values.get_function_values(rel_solution, solution_values);
   fe_values.get_function_gradients(rel_solution, solution_grads);

int cell_index = cell->active_cell_index();

//  std::cout << cell_index << std::endl;

Tensor<2, dim> grad_u, eps;

for (unsigned int q = 0; q < n_q_points; ++q)
   {
 std::cout << "\n - in energy GAUSS point " << q + 1 << " - 
\n" << std::endl;

  for (unsigned int i = 0; i < dim; ++i)
 {
std::cout << "GRADIENT 1: " << solution_grads[q][0] 
<< std::endl;
std::cout << "GRADIENT 2: " << solution_grads[q][1] 
<< std::endl;

 // Compute displacement gradients
grad_u[i] = solution_grads[q][i];

 // Compute strains
 eps = 0.5 * (grad_u + transpose(grad_u));

  }

  print_tensor(grad_u, "DISPLACEMENT GRADIENT");
  }

  }
}

The above code gives me the following output:

- in GAUSS point 1 - 

GRADIENT 1: -0.000312
GRADIENT 2: -0.001801
GRADIENT 1: -0.000312
GRADIENT 2: -0.001801

DISPLACEMENT GRADIENT
 0.001081   -0.001477  
 -0.00   -0.008525  


STRAIN TENSOR
 0.001081   -0.000738  
 -0.000738   -0.008525  


 - in GAUSS point 2 - 

GRADIENT 1: 0.000312
GRADIENT 2: -0.001801
GRADIENT 1: 0.000312
GRADIENT 2: -0.001801

DISPLACEMENT GRADIENT
 0.001081   0.001477  
 -0.00   -0.008525  


STRAIN TENSOR
 0.001081   0.000738  
 0.000738   -0.008525  


 - in GAUSS point 3 - 

GRADIENT 1: -0.001164
GRADIENT 2: -0.006723
GRADIENT 1: -0.001164
GRADIENT 2: -0.006723

DISPLACEMENT GRADIENT
 0.004034   -0.001477  
 -0.00   -0.008525  


STRAIN TENSOR
 0.004034   -0.000738  
 -0.000738   -0.008525  


 - in GAUSS point 4 - 

GRADIENT 1: 0.001164
GRADIENT 2: -0.006723
GRADIENT 1: 0.001164
GRADIENT 2: -0.006723

DISPLACEMENT GRADIENT
 0.004034   0.001477  
 -0.00   -0.008525  


STRAIN TENSOR
 0.004034   0.000738  
 0.000738   -0.008525  


And the below are results from the postprocessor computation for a single 
element (already checked them, they should be correct):

- in GAUSS point 1 - 


DISPLACEMENT GRADIENT
 0.00   -0.002557  
 0.00   -0.008525  


STRAIN TENSOR
 0.00   -0.001279  
 -0.001279   -0.008525  


 - in GAUSS point 2 - 


DISPLACEMENT GRADIENT
 0.00   0.002557  
 0.00   -0.008525  


STRAIN TENSOR
 0.00   0.001279  
 0.001279   -0.008525  


 - in GAUSS point 3 - 


DISPLACEMENT GRADIENT
 0.005115   -0.002557  
 -0.00   -0.008525  


STRAIN TENSOR
 0.005115   -0.001279  
 -0.001279   -0.008525  


 - in GAUSS point 4 - 


DISPLACEMENT GRADIENT
 0.005115   0.002557  
 -0.00   -0.008525  


STRAIN TENSOR
 0.005115   0.001279  
 0.001279   -0.008525  

Why is the fe_values.get_function_gradients() function giving me wrong 
gradients?

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: Relationships between material parameters

['Seyed Ali Mohseni' via deal.II User Group]

Your ideas are very constructive and useful, I didn't think it would be so 
complex to make it general in deal.II.
I am not an expert still, but how about we make it something similar to the 
FEValues class where you just specify the necessary parameters and the 
desired conversion.

@ Jean-Paul: I am now a bit busy until February, but it would be a great 
honor for me to contribute to deal.II, if possible.

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Relationships between material parameters

['Seyed Ali Mohseni' via deal.II User Group]

Hi everyone,

I am just implementing some material parameter relationships such as  

lambda = E * v / (1 + v) / (1 - 2 * v)

for instance.

There is a table for such conversion of quantities.

Would it be beneficial to have such a library that computes and converts 
these quantities automatically in Deal.II?
Or is there already such a library that I can access instead of writing 
everything myself?

If there is not, I would be willing to implement such a library and add it 
to the core of Deal.II, so people can access it.

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.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Relationship between material parameters

['Seyed Ali Mohseni' via deal.II User Group]

Hi everyone,

I am just implementing some material parameter relationships such as  

lambda = E * v / (1 + v) / (1 - 2 * v)

for instance.

There is a table for such conversion of quantities.

Would it be beneficial to have such a library that computes and converts 
these quantities automatically in Deal.II?
Or is there already such a library that I can access instead of writing 
everything myself?

If there is not, I would be willing to implement such a library and add it 
to the core of Deal.II, so people can access it.

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.
For more options, visit https://groups.google.com/d/optout.