Steps without NLGEOM that ccx iterates (e.g. an unreached *PLASTIC): *RIGID BODY, *KINEMATIC and mean-rotation MPCs use finite rotations

Summary

Section 6.7.4 (Kinematic and Distributing Coupling) of the 2.23 manual says:

Coupling constraints generally lead to nonlinear equations. In linear calculations (without the parameter NLGEOM on the *STEP card) these equations are linearized once and solved. In nonlinear calculations, iterations are performed in each of which the equations are linearized at the momentary solution point until convergence.

ccx also solves a step without NLGEOM with Newton-Raphson iterations when the model has a nonlinear material such as *PLASTIC, or a contact pair. The output then says “Newton-Raphson iterative procedure is active”, but not “Nonlinear geometric effects are taken into account”, so the elements stay geometrically linear. In such a step ccx re-evaluates the nonlinear MPCs (*RIGID BODY, *KINEMATIC coupling, the mean-rotation MPCs of shell and beam rotations) at the current displacements in every iteration. They therefore describe finite rotations, as with NLGEOM.

The result is that a *PLASTIC card whose yield stress is never reached changes the displacements. Test: a rigid face turned 0.3 rad about z.

deck node 11 (x, y)
*ELASTIC only -0.1500000, 0.1500000 (linear: theta x r)
the same plus *PLASTIC with yield stress 1e9 -0.1700919, 0.1254283 (finite: (R - I) r)
*ELASTIC only, *STEP, NLGEOM -0.1700919, 0.1254283

The material never yields, so the first two decks describe the same problem and should give the same answer. They differ by 13 % in x and 16 % in y. A *KINEMATIC coupling of the same face gives the same three results. Exit 0 in all cases; the only *WARNING is the default minimum-increment one, which the elastic deck prints too.

Where the elements carry the load, the answer of the iterated step is neither the linear nor the NLGEOM one. An S4 strip under a tip moment (a mean-rotation MPC) that turns its tip by 0.2 rad deflects 4.999500E-02 with *ELASTIC (beam theory: 0.049995), 4.982203E-02 with NLGEOM, and 4.818387E-02 with the unreached *PLASTIC (the deck is below).

Version

ccx 2.23 (macOS arm64) and current master (537aa81), same numbers.

Already reported or fixed?

Not reported: I found no matching GitHub issue (open or closed) or pull request in Dhondtguido/CalculiX, and no calculix.discourse.group topic (searched 3 Nov 2026). Not fixed: the nonlinmpc calls in nonlingeo.c are unchanged on master (537aa81).

Minimal input

rb_plastic.inp:

** Two C3D8 in a column (z = 0 to 2), bottom face held. The top-face nodes 9-12
** form a *RIGID BODY (REF NODE 100, ROT NODE 101 at the face centre), which
** turns 0.3 rad about z. No NLGEOM. The material has a *PLASTIC card whose
** yield stress (1e9) is never reached.
*NODE, NSET=NALL
1, 0., 0., 0.
2, 1., 0., 0.
3, 1., 1., 0.
4, 0., 1., 0.
5, 0., 0., 1.
6, 1., 0., 1.
7, 1., 1., 1.
8, 0., 1., 1.
9, 0., 0., 2.
10, 1., 0., 2.
11, 1., 1., 2.
12, 0., 1., 2.
100, 0.5, 0.5, 2.
101, 0.5, 0.5, 2.
*ELEMENT, TYPE=C3D8, ELSET=EALL
1, 1, 2, 3, 4, 5, 6, 7, 8
2, 5, 6, 7, 8, 9, 10, 11, 12
*NSET, NSET=BOT
1, 2, 3, 4
*NSET, NSET=TOP
9, 10, 11, 12
*MATERIAL, NAME=M
*ELASTIC
1000., 0.3
*PLASTIC
1.e9, 0.
*SOLID SECTION, ELSET=EALL, MATERIAL=M
*RIGID BODY, NSET=TOP, REF NODE=100, ROT NODE=101
*BOUNDARY
BOT, 1, 3
100, 1, 3
101, 1, 2
101, 3, 3, 0.3
*STEP
*STATIC
0.1, 1.
*NODE PRINT, NSET=TOP
U
*END STEP

Variants:

  • rb_elastic.inp: the two lines *PLASTIC / 1.e9, 0. deleted.
  • rb_nlgeom.inp: rb_elastic.inp with *STEP, NLGEOM instead of *STEP.
  • kin_plastic.inp / kin_elastic.inp: the same with a *KINEMATIC coupling instead of the *RIGID BODY:
-101, 0.5, 0.5, 2.
 ...
-*RIGID BODY, NSET=TOP, REF NODE=100, ROT NODE=101
+*SURFACE, NAME=STOP
+2, S2
+*COUPLING, REF NODE=100, SURFACE=STOP, CONSTRAINT NAME=C1
+*KINEMATIC
+1, 3
 *BOUNDARY
 BOT, 1, 3
-100, 1, 3
-101, 1, 2
-101, 3, 3, 0.3
+100, 1, 5
+100, 6, 6, 0.3
strip_plastic.inp: S4 strip with a tip moment (mean-rotation MPC)
** S4 strip 0.5 x 0.1 x 0.01 clamped at x = 0 (nodes 1, 2); a moment about -y
** is shared by the tip nodes 11, 12. The material has a *PLASTIC yield stress of 1e9 (never reached).
*NODE, NSET=NALL
1,0,0,0
2,0,0.1,0
3,0.1,0,0
4,0.1,0.1,0
5,0.2,0,0
6,0.2,0.1,0
7,0.3,0,0
8,0.3,0.1,0
9,0.4,0,0
10,0.4,0.1,0
11,0.5,0,0
12,0.5,0.1,0
*ELEMENT,TYPE=S4,ELSET=P
1,1,3,4,2
2,3,5,6,4
3,5,7,8,6
4,7,9,10,8
5,9,11,12,10
*MATERIAL,NAME=M
*ELASTIC
2e5,0.
*PLASTIC
1.e9, 0.
*SHELL SECTION,ELSET=P,MATERIAL=M
0.01
*NSET,NSET=TIP
11,12
*BOUNDARY
1,1,6
2,1,6
*STEP
*STATIC
0.1, 1.
*CLOAD
11,5,-0.0003333
12,5,-0.0003333
*NODE PRINT,NSET=TIP
U
*END STEP

Node 11 at the end of the step: (1.0E-17, 7.6E-17, 4.818387E-02). Without the *PLASTIC lines: (8.4E-16, 5.4E-15, 4.999500E-02); beam theory M L² / (2 E I) = 0.049995. Without them and with *STEP, NLGEOM: (-3.326140E-03, ~0, 4.982203E-02); the elastica gives (-3.3266E-03, 0.049832).

How to run

ccx rb_plastic, ccx rb_elastic, ccx rb_nlgeom (and the kin_ variants); compare node 11 in the .dat files.

Expected

Node 11 at (-0.1500000, 0.1500000) in both steps without NLGEOM, as the manual sentence above describes (“linearized once”), the same as rb_elastic. The constraint should not depend on whether the material model has a *PLASTIC card that never becomes active.

Actual

deck node 11 (x, y, z) solver output
rb_elastic, kin_elastic -1.500000E-01, 1.500000E-01, 0 linear
rb_plastic, kin_plastic -1.700919E-01, 1.254283E-01, 0 “Newton-Raphson iterative procedure is active”, no “Nonlinear geometric effects” line
rb_nlgeom -1.700919E-01, 1.254283E-01, 0 “Nonlinear geometric effects are taken into account”

Exact rigid rotation of r = (0.5, 0.5, 0) by 0.3 rad: (R - I) r = (-0.1700919, 0.1254283, 0); linear: theta x r = (-0.15, 0.15, 0).

Where it goes wrong

linstatic.c never calls nonlinmpc, so in a linear step the nonlinear MPCs keep the linearisation made when they were generated (at the undeformed state). nonlingeo.c calls nonlinmpc with the current displacements vold at the start of the step (line 951), in every increment (line 2071) and in every iteration (line 2442), whether or not the step is geometrically nonlinear (iperturb[1]).

ccx’s own definition of the flag: CalculiXstep.c:134, “iperturb(2) (ITG): if 0: geometrically linear calculation (linearized strain)”. And for contact, contactpairs.f:244-247: “the definition of a contact pair triggers a call to nonlingeo (for static calculations) but not automatically to the nonlinear calculation of strains (i.e. iperturb(2) should be zero unless NLGEOM is activated)”. The MPC update ignores that distinction. (Section 10.7.2 (Nonlinear calculations) lists “update the coefficients of nonlinear MPC’s, if any” in every iteration of nonlingeo, which describes what the code does; it does not say whether that is meant for steps without NLGEOM, and for coupling constraints section 6.7.4 says it is not.)

Fix

When iperturb[1] is 0, keep the coefficients of the nonlinear MPCs at the undeformed state, as linstatic.c does, e.g. by evaluating them once and skipping the later nonlinmpc updates in nonlingeo.c for such a step. (Not tried.) If iterating them is intended, section 6.7.4 should say that the coupling equations are iterated whenever the Newton-Raphson procedure is active, not only with NLGEOM.

Confirmed using official ccx 2.23 build on Windows, needs a GitHub ticket.

rb_plastic.dat:

                                INCREMENT     6


 displacements (vx,vy,vz) for set TOP and time  0.1000000E+01

         9  1.700919E-01 -1.254283E-01  0.000000E+00
        10  1.254283E-01  1.700919E-01  0.000000E+00
        11 -1.700919E-01  1.254283E-01  0.000000E+00
        12 -1.254283E-01 -1.700919E-01  0.000000E+00