Global constraints via Lagrangian Multiplyers [feap v8.6]

[a118145]

Dear feap forum,

I am currently experimenting with the implementation of global constraints via Lagr. Mult., see programmer manual 8.6 section 5.11. Following the instructions there (which is more or less identical to the previous feap programming manual 8.5), I do not achieve a working solution. However, I stumbled upon additional comments in section 5.5.1 regarding global constraints. If I solely set

Code: [Select]

GLOBal
    EQUAtion,,3  ! three global constraint equations

in the input file (directly after first FEAP line) as implied by sec. 5.11, the implementation does not work: The parameter nge, which is mentioned in sec. 5.5.1 is still zero, but the other integer geqnum, which indicates the number of global equations according to sec. 5.11, and hence nst are set properly. Due to nge being zero, the corresponding solution array ulg is not updated, either. It seems to me as if nge is mandatory. Indeed, if i set nge=3 during task 1 in my user element, it works fine - at least with the internal feap solver. Doing so, I can even omit the GLOBal statement in the input file at all and only rely on setting nge=3 during task 1, and still get a working solution.

I am quite confused how to set the necessary integer values correctly. Can you clarify this?

Warm regards,

a118145

[Prof. R.L. Taylor]

1. In your user element you need to include 'hdata.h' -- which has the parameter 'nge'
2. In each element input the value you use for 'nge' as part of isw.eq.1
3. During isw.eq.3 (and isw.eq.6) you need to compute the residual and tangent for the global equations (see pp 53 of Programer Manual for how element array is organized). Check that "nst" has been adjusted to be = ndf*nen + geqnum and add your development to last row/cols.

You can check if the global equations are set looking at the array "IE(5,*)" where * is the material set.

You can check if the solver knows there are global equations using "SHOW" and "SHOW ID" - show should have an output for number of equations and "ID" has one part that contains the equation numbers for each node.  The last node should have a value that is NEQ - GEQNUM.



You may not need to set the global command if all elements have contribute to all global equations. 

[Prof. R.L. Taylor]

You will need the global command, the parameter 'geqnum' i used to update the solutions in ulg(*), so without it you do not get an update.f

[a118145]

Dear Prof. Taylor,

Thank you for your clarification! I cross-checked the remarks 1-3 and the status of the global equations in IE(5,*). Everything seems fine. Regarding your second post about the necessity of setting the GLOBal option in the input file, I am not sure about this. I tried two identical setups but commented out the option in one case: ulg seems to be updated either way.

Regarding the SHOW and SHOW ID options, I observe the following. Note, that I switched to a two-dimensional 4-node unit square single element test. As I have two displacements at each node, the SHOW output seems correct to me:

Code: [Select]

     Mesh Dimension    =           2  :  Number Dof/Node   =           2
     Number Nodes      =           4  :  Number Elements   =           1
     Number Materials  =           1  :  Number Equations  =           8
     No. Global Eqs.   =           2  :  Global Partition  =           0

The number of global equations is set correctly, too. I find it somewhat difficult to interpret the output of SHOW ID:

Code: [Select]

--> Part  1 of array

    Matrix: ID                           
    row/col         1         2         3         4
       1            1        -1        -2         5
       2            2         3         4         6
 --> Part  2 of array

    Matrix: ID                           
    row/col         1         2         3         4
       1            0         1         1         0
       2            0         0         0         0

How do I interpret the eight numbers in Part 1 and Part 2, which are somehow linking the node number (first row) and the glob. eq. number (first column)? I suppose, they indicate a position or position offset in the global tangent matrix. The global tangent is 10x10 in my case (single 4-node 2D element + 2 glob. eq. --> cross-checked via nst) and the according submatrices are set in the user element according to Figure 5.7 in the most recent programmer manual.

Interestingly, I get an error in the routine INZA as soon as I try to employ the Pardiso solver and not the feap internal one:

Code: [Select]

    INZA: Eq. No. =       7 Col. No. =       1
          Col. top =      22 Col. Bot =      -1
          Column entries

 --> ERRORS OCCURRED: For details see file: O_1elem001

Is this error connected to a faulty indexing, that maybe becomes apparent from my SHOW ID output?

Thank you very much for your continuous effort and support!

Warm regards,

a118145

[Prof. R.L. Taylor]

1. The value of 'geqnum' will be set form the element values of 'nge' if not input by the global command.

2. The Part 1 of ID are the equation numbers.  Negative or zero values are ignored as they indicate the dof will be set by a boundary condition.

3. Part 2 of the ID array are boundary condition values.  A zero indicates an active equation location and a non-zero indicates a specified b.c. dof.  In many problems the data is a bit redundant, however for problems with partitions it is used independently.

4. The error using pardiso may be a bug, probably related to global equations -- needs to be checked more.

[Prof. R.L. Taylor]

I have confirmed that neither using Pardiso or the feap sparse solver works for problems with global equations.

[a118145]

Dear Prof. Taylor,

Thank you for figuring this out. This is not too good for me. Let me ask you three more questions:

1) Based on your observations, what has to be done in order to make this work for the Pardiso solver?

2) Are the global constraints still not available in parfeap? Actually, this is essential for me. If not, what would be the starting point of enhancing the feap 8.6 implementation?

3) Is there a way to prescribe the Lagrange Multiplyers instead of the RHS - similar to specifying Dirichlet instead of Neumann BCs?

Warm regards,

a118145

[Prof. S. Govindjee]

We are working on a patch to allow this.

[epn]

Dear Prof. Govindjee, dear Prof. Taylor,

I am also currently trying to implement a global constraint. The above discussions are quite helpful.

I would like to ask if the problem with the sparse solvers is resolved?

Thank you for the information in advance!

Best regards