Compute Compliance w/ Prescribed Displacements

[jbruss]

Hello -

I have a prescribed displacement field and no external forces/tractions. Is there a way for me to compute the compliance in this case? Before when I was applying a force load, I could use F (hr(np(27))) and U (hr(np(40))) to compute this quantity (F^T U) but now F does not contain the reaction information at the dirichlet degrees of freedom. Is this information somewhere in hr(np(26)) after a solve or do I need to put the solution, U, into a PETSc Vec and perform the matrix scalar product with the current tangent matrix (i.e. U^T K U)?

Thank you for your help!
Jonathan

[Prof. S. Govindjee]

Did you try using REACtion?  It should compute the reactions at the driven degrees of freedom.

[jbruss]

Hi Professor -

If I call REACtion, are the reactions placed in the first ndf*numpn entries of hr(np(26))? I see that this vector is changed but I don't know how to interpret the quantities within it.

Thank you again for your help.
Jonathan

[FEAP_Admin]

To learn what built in commands do you can follow these steps.

1) Look in program/pmacr.f at the data statement for 'wd'.

2) Find the command of interest in the data statement.  The number in the continuation column will will tell you which pmacrX.f file contains the code for the command.

3) Search now in the identified pmacrX.f file for your command.  Easiest is to search for '[cmd-name', where cmd-name is the first 4 letters of the command of interest.

[jbruss]

Hi -

Thank you for your help. I see now what information is in this vector but it does not contain the information I am looking for. When the degrees of freedom in the system with Dirichlet constraints are condensed out of the system, using the command REACtion computes the nodal reactions at every node except for these degrees of freedom (for which it reports a value of 0.0, regardless of whether the prescribed displacement is nonzero).

I confirmed this for an extremely simple problem with 2 quadrilateral elements in 2D. The left side is fixed and the right side I prescribed a displacement of 0.001 (the length of the "bar" is 2.0 and the width is 1.0). I obtain the correct displacement field. However using REACtion, the reaction values at the left and right ends are given as 0.0 and I get nonzero values for the degrees of freedom that are not constrained.

Is there someway to obtain the reactions for this problem at nodes which are assigned a prescribed value if these degrees of freedom are condensed out of the system?

Thank you again,
Jonathan

[Prof. R.L. Taylor]

Can you tell us if you are using the parallel code or the serial code when you do the reaction which has the error?

[Prof. S. Govindjee]

Also, if you are doing this in parallel, which OUTDomain scheme are you using?  And which version of the program are you using?

[jbruss]

Hi Professor Taylor and Govindjee -

Unfortunately, I am using the parallel version of FEAP 8.3. I'm sort of stuck with using 8.3 since the code is heavily patched and I'm really not sure how much effort it would take to migrate all of the existing code our group and I have written into 8.5. The outdomain scheme is the default AIJ and I am executing the code on only one processor. Perhaps this error is fixed in FEAP 8.5.

In order to find a way around this, I simply created a new isw in the element I am using and sum the internal force (since in quasi-static simulations this should equal the external force). This works well enough for me for now. Eventually I may try to migrate the code to 8.5 because I would really like to be using the newest version of FEAP (and of PETSc). Unfortunately, I don't have the time at the moment.

Are the reactions computed correctly at Dirichlet boundaries in 8.5? I'm sure there have been many enhancements and bug fixes since 8.3 so I will do my best in the next few months to attempt an upgrade.

Thank you again for your time,
Jonathan

[FEAP_Admin]

I checked version 8.5 in parallel on a simple tension strip, clamped at x=0 and having an applied force at x=L.  The command reac,coor,1,0 produced the correct reaction forces at the clamped end.

I can also see in the comments log in parfeap/pform.f that this problem was apparently corrected on 13/12/2013, but I can not see exactly how to back port the change to version 8.3.

You may wish to consider upgrading to the latest version.

[jbruss]

Perfect. That is what I needed to know! Thank you all!

[Prof. R.L. Taylor]

One other comment.  The way you proposed to compute the reactions with a new isw may not work when you have more domains since we also include additional "ghost" elements in each domain to minimize data exchange.

[jbruss]

Professor Taylor -

Once I get things working in serial I will try to modify it to work in parallel and check for ghost elements. Thank you for the heads up.

Jonathan