This example shows a rudimentary script for the simulation of contact. It can definitely be enhanced and we are open to specific improvements.

Here, we're using a beam that is clamped on one side. The beam is largely rigid, but near the clamping, it is flexible. The flexible part acts like a hinge. We're using an actuation force of 1 N in x-direction on node 3. Without any contact, the tip of the beam would displace 0.0311 m.

We're imagining a line at x=0.025, representing a contact surface that cannot be crossed. To simulate this, we're creating a custom force on the node of contact (node number 4), whenever it crosses the position of x=0.025. We choose the value of this force to be -stiffn*(pos-x0). This creates a fictitious force that acts like a spring with stiffness stiffn and reference position of x0=0.025. Variable pos is the current value of the suspected node of contact.

This simulation is tricky because it involves iterations to get the system in the right state. We currently have no clear rationale for picking the values of stiffn of the fictitious force (or the precise form of this force, for that matter). Also, it's best to do these iterations within one single Spacar call, seeing as each new SPACAR Light call will make the system start from the undeformed configuration again. The current interface for doing a custom force within Spacar is by means of the usersig interface. This interface hasn't been implemented in SPACAR Light properly, but here are some pointers to get a simulation up and running.

• Get the development version of SPACAR Light that supports the opt.customdyn field
• Use opt.customdyn = 'usersig PUSHSIG' to have Spacar call the m-file PUSHSIG.m (all of this is case-sensitive)
• This m-file will be called repeatedly, depending on the number of opt.loadsteps. This file contains a function with very specific input and output requirements. Within this function, you'll have available the x-vector, which contains the current positions of the nodes. Unfortunately, the node numbering in Spacar (full) is different from SPACAR Light. This is actually by design, and something you usually shouldn't have to come across. Here, you do, though. The conversion is simple and detailed in the file. Note that this simple conversion may not work whenever you use multiple elements per beam with eprops(i).nbeams > 1. (The .dat-file actually shows you what the conversion should be.)
• debug
• loadsteps ignored if no force
• global

PUSHSIG.m

function [time,sig,f] = PUSHSIG(t,ne,le,e,ep,nx,lnp,x,xp)
 
    x0 = 0.025; %defines contact line
    stiffn = 300;
    spalightnode = 4;
 
    %conversion between SPACAR Light node and Spacar (full) node.
    %(for each SPACAR Light node, there are 2 nodes in Spacar full.
    %They are numbered consecutively. First one is the translation node,
    %second one is the rotation node. Here, we're taking just the translation node.)
    spanode = (spalightnode-1)*2+1;
 
    %get current position
    pos = x(lnp(spanode,1:3)); %contains x,y,z coordinate of spalightnode
 
    %define a custom force
    node = spanode; %node number (according to Spacar full numbering)
    dir = 1; %direction of force: 1 for x, 2 for y, 3 for z
 
    debug_info = 1; %print out some debug info during simulation
 
    if debug_info
        fprintf('Current position: %f\n',pos(1))
    end
 
    if pos(1) > x0 %if within (forbidden) contact region
        val = -stiffn*(pos(1)-x0); %define force to push it back
        deriv = stiffn;
        f = [node dir val deriv 0]; %apply force
 
        if debug_info
            fprintf('Within contact region. Applying force: %f\n',val)
        end
 
    else
        f = []; %if not within contact region, don't apply contact force
        if debug_info
            fprintf('Not within contact region.\n')
        end
    end
 
    sig = []; %can be left empty (is used for specifying generalized stresses)
    time = t; %this line needs to be there; don't modify.
 
end