How to output a single scalar in a parallel code

[Hamed Babaei]

Hello,

I want to make an output of the resultant force existing on external surfaces of a cube domain in a parallel code. what I did in the serial code was that I integrated stress on external faces and as a result, I got one scalar resultant force for the face to be written in a text file in every time step. However in parallel, a face may be possessed by several processors each of which producing a portion of the overall force existing on the whole face. I was wondering what may be the best approach in parallel codes to make an output file of one single scalar output which is obtained by integration over domain. 

Best regards,

Hamed

[RAJAT ARORA]

Hello,

I am solving a similar problem. What I do is to sum the contributions from all the processors, and then use Utilities::MPI::Sum() to get the overall force which I then write in a file on a master process,

For the contribution of a single cell, you have to integrate stress X normal X area for each locally owned cell on faces on the boundary.

[Hamed Babaei]

Dear Rajat,

 

... which I then write in a file on a master process.

I was wondering how to write in a file on a master process so that there would be just one output file not as many as processor exist. I use the following commands to write in my file in serial code :

std::ofstream myfile;

myfile.open ("resultant_strass.txt");

myfile<< resultant_stress<<std::endl;

myfile.close();

Do I need to make any change in it?

Thanks

[Timo Heister]

> Do I need to make any change in it?

You need to only do it on processor with rank 0, for example how we do
it in step-40:
https://urldefense.proofpoint.com/v2/url?u=https-3A__github.com_dealii_dealii_blob_master_examples_step-2D40_step-2D40.cc-23L610&d=DwIBaQ&c=Ngd-ta5yRYsqeUsEDgxhcqsYYY1Xs5ogLxWPA_2Wlc4&r=4k7iKXbjGC8LfYxVJJXiaYVu6FRWmEjX38S7JmlS9Vw&m=lQtnDi-x-kIhyIM5zjnGN6iYg0vOkjsqRYUpdu5AqQ0&s=gZQBXRI7J6MIAvVrgIt6OQ9GPhPQXm2yMlR6AnwcNVE&e=

--
Timo Heister
http://www.math.clemson.edu/~heister/

[RAJAT ARORA]

Hello,

Please look at the solution by Timo. 

In your case, you can do something like.

 if (Utilities::MPI::this_mpi_process(mpi_communicator) == 0)
{

[amir.pa...@gmail.com]

Dear Hamed,

Integration stress on external faces is correct when the material constitutive model is linear (for example elastic material model in which behavior is not history dependent). It is better to extract residual forces at the boundary from the system residual vector. 

Regards

Pasha

[Jean-Paul Pelteret]

Dear Amir (and Hamed),

 

Integration stress on external faces is correct when the material constitutive model is linear (for example elastic material model in which behavior is not history dependent). It is better to extract residual forces at the boundary from the system residual vector. 

I disagree and consider both statements that you've made to be incorrect. For starters, in a nonlinear problem one uses a nonlinear solution scheme to drive the residual (a consequence of equilibrium imbalance) down to zero. So not only does the residual not represent any measure of external force, but also one expects this vector to be zero in for a balanced system!

In general (for a linear/nonlinear hyperelastic/dissipative material), to compute the total force acting on the surface of a domain one should integrate the tractions acting on that surface. As the body is in equilibrium, they are balanced by the stresses generated by the material. Using Cauchy's theorum one can re-express this traction in terms of the internal stresses and compute the forces developed on any arbitrary cut-plane in the material domain.

To demonstrate, I've attached an example for one-field elasticity that shows the computation of the reaction forces on either side of a cube under uniaxial stress conditions. All of the computations are done in the post-processing step (Solid<dim>::output_results). Since this cube has flat faces, computing the nodal reaction forces (which, by the way, are outputted so that you can visualise them as well) and then summing them is equivalent to directly integrating the contributions over the surface.

I hope that this helps provide you with some clarification on the point that you made.

Regards,

Jean-Paul

[amir.pa...@gmail.com]

Dear Jean-Paul

Why we must calculate reaction forces at a boundary which has external force? When we apply external load (load control solution) so we control the load level at that boundaries and reaction forces must be equal to applied load level (residual is zero) at each time step at these BCs due to equilibrium condition. But we usually calculate reaction forces at BCs with prescribed displacement which it is better to calculate them from residual vector (which equal to [-1.0 * internal forces] at these BCs). You aslo explained same method as:

"In general (for a linear/nonlinear hyperelastic/dissipative material), to compute the total force acting on the surface of a domain one should integrate the tractions acting on that surface. As the body is in equilibrium, they are balanced by the stresses generated by the material."

As the body is in equilibrium, they are balanced by the stresses generated by the material = internal forces

So when we calculate internal forces why computation is performed again for integral over BCs? 

I think this method is used in advance FEM software. 

Please, did you test the attached code for a plastic material model?

Best Regards

Pasha

[Jean-Paul Pelteret]

Dear Amir,

 

Why we must calculate reaction forces at a boundary which has external force?

You certainly don't have to. In my example I used this as a demonstration that the forces computed from the internal stresses matched the applied (uniformly distributed) load. If you believe this then you could happily apply the same procedure to any surface undergoing displacement control, knowing that it will produce the desired result.

 

But we usually calculate reaction forces at BCs with prescribed displacement which it is better to calculate them from residual vector (which equal to [-1.0 * internal forces] at these BCs). 

As the body is in equilibrium, they are balanced by the stresses generated by the material = internal forces

Ok, I see now where you're coming from. This might be true, but only for constrained degrees-of-freedom for quasi-static problems. But in his initial post Hamed wasn't specifically referring to reaction forces (unless I misunderstood something), but rather the resultant force on an unspecified boundary surface. Furthermore, for time-dependent problems you may have additional contributions to these components of the residual due to the imposition of constraints on the velocity and acceleration field (if these constraints are non-homogeneous).

 

Please, did you test the attached code for a plastic material model?

No I haven't. Do you anticipate the result will be any different when using a plasticity model? If so, how?

Regards,

Jean-Paul

[amir.pa...@gmail.com]

Dear jean-Paul

Yes, we can use internal forces for all situations. I must correct my answer by replacing "residual vector" by "internal force vector". However, can we calculate internal force by integrating over boundary surfaces when material model is history dependent and has internal state variable such as plasticity? I think the answer is no, so I have suggested to use "internal force vector" in the first answer instead of integration over the boundary surfaces.

Thank you,

Pasha

[Jean-Paul Pelteret]

Dear Pasha,

 I must correct my answer by replacing "residual vector" by "internal force vector"

Ok, great. So this is a slightly different situation. I was going to say that I think that we're cross-talking here, and you saved me from sketching out a whole bunch of equations to clarify things :-) 

Yes, we can use internal forces for all situations. 

Hmm... I still don't think so. What about if you have body forces? Then the internal stresses balance both the externally applied loads plus the volume forces. Do you agree with this? However, one can of course compute the internal stress plus the body force loading... 

But then consider a cube with a non-uniform pressure on more than one side. How do I compute from this global vector the average or total load applied to any individual surface? You can't, because the values in the degrees-of-freedom shared by any two sides of the geometry contain load contributions from two faces. For reaction forces at constrained surfaces - no problem of course!

Also, I think that one should consider which vectors one typically constructs during assembly. The most efficient way to do assembly (but I recognise that this is not universally followed - not even by me) is to compute the constrained RHS vector and tangent matrix in one shot during assembly. That means that you never need to actually construct the individual internal, external and body force vectors in order to solve your linear system because, well, they get bundled into the RHS and then further manipulated when imposing constraints.

Dependent on what it is you want to achieve, I do concur that what you suggest may be perfectly plausible. But in general... I'm still unconvinced. To be introspective, what I suggest (and what Hamed has done here) might be overkill for plain elasticity alone since you can always recompute the traction vector itself as it is you that specified the loading at each time-/load-step. In coupled problems, however, external tractions are not always so easily defined and its actually sometimes easier to simply do this stress-type integration procedure.

However, can we calculate internal force by integrating over boundary surfaces when material model is history dependent and has internal state variable such as plasticity

Yes, you can - its just not very convenient. You could do this in one of two ways:

1. Extrapolate the internal variables to the faces and recompute the stresses with these values. Its clear that this might be problematic because one has to make an assumption about whether there is plastic yield at the extrapolated points (as can be attested to by the note in the ContinuousQuadratureDataTransfer class - we'd like to extend this at some point and we've already worked out a couple of options that may work for plasticity but simply haven't tried to implement them yet), so I would recommend option 2...

2. Simply create and track additional quadrature point data (with internal data) at the faces at which you'll be doing this post-processing. That way you can update the internal variables as necessary and they'll always be in sync with the solution field. You can then compute stresses etc. at these points whenever necessary. Its annoying if the local problem is expensive to compute, but the number of face QP data that you'd create is likely far less than those required for assembly etc.

Best,

Jean-Paul