Simufact Forming 2023.

2
Scientific Tutorial

User subroutines, User-defined kinematic,

User-defined material
Copyright © 2023 Hexagon AB and/or its subsidiaries. All rights reserved.
This documentation, as well as the software described in it, is furnished under license and may be used only in accordance with the terms of such
license.

Hexagon reserves the right to make changes in specifications and other information contained in this document without prior notice. The concepts,
methods, and examples presented in this text are for illustrative and educational purposes only, and are not intended to be exhaustive or to apply
to any particular engineering problem or design. Hexagon assumes no liability or responsibility to any person or company for direct or indirect
damages resulting from the use of any information contained herein. This notice shall be marked on any reproduction of this documentation, in
whole or in part. Any reproduction or distribution of this document, in whole or in part, without the prior written consent of Hexagon is prohibited.

This software may contain certain third-party software that is protected by copyright and licensed from third party licensors. Additional terms
and conditions and/or notices may apply for certain third-party software. Such additional third-party software terms and conditions and/or notices
may be set forth in documentation and/or at http://www.mscsoftware.com/thirdpartysoftware (or successor website designated by Hexagon from
time to time).

This product includes software developed by the University of Chicago, as Operator of Argonne National Laboratory.

PCGLSS 8.0, Copyright © 1992-2023 Computational Applications and System Integration Inc. All rights reserved. Portion of this software are
owned by Siemens Product Lifecycle Management Software, Inc. © Copyright 2023.

Trademarks

The Hexagon logo, Hexagon, MSC, MSC Nastran, Adams, Dytran, Marc, Mentat, Patran, Simufact, Simufact Forming, Simufact Welding and
Simufact Additive are trademarks or registered trademarks of Hexagon AB and/or its subsidiaries in the United States and/or other countries.

Windows is a trademark or registered trademark of Microsoft Corporation in the United States and other countries.

Linux is a registered trademark of Linus Torvalds.

3Dconnexion and SpaceDevice are registered trademarks of 3Dconnexion Inc.

Nastran is a registered trademark of NASA.

FLEXlm and FLEXNet Publisher are trademarks or registered trademarks of Flexera Software.

Parasolid is a registered trademark of Siemens Product Lifecycle Management Software, Inc.

Magics and Materialize are trademarks or registered trademarks of Materialize NV and/or its subsidiaries.

All other trademarks belong to their respective owners.

Contact:

http://www.simufact.com/contact
2023.2 Simufact Forming

Table of Contents
1. User subroutines ........................................................................................................ 1
1.1. FORTRAN guide (for FE-Solver) ....................................................................... 2
1.1.1. Common block ..................................................................................... 3
1.1.2. File operations ...................................................................................... 5
1.1.3. Conditional statement ............................................................................ 6
1.1.4. Loops .................................................................................................. 7
1.1.5. Domain decomposition method (DDM) ..................................................... 8
1.1.6. Debugging ......................................................................................... 10
1.2. Variable units ................................................................................................ 12
1.3. Coordinate axis ............................................................................................. 13
1.4. User subroutine and utility routines ................................................................... 13
1.4.1. Subroutine details ................................................................................ 14
1.4.2. Subroutine integration .......................................................................... 18
1.4.3. User variables, state variables and post variables ....................................... 19
1.4.4. Utility routines .................................................................................... 21
1.4.5. Common user subroutines ..................................................................... 28
1.4.6. Post codes .......................................................................................... 30
1.5. User subroutine: Examples .............................................................................. 32
1.5.1. State and post variables ........................................................................ 32
1.5.2. Motion .............................................................................................. 41
1.5.3. Damage ............................................................................................. 45
1.5.4. Friction .............................................................................................. 51
1.5.5. Table evaluation .................................................................................. 59
1.5.6. Restart with global variables ................................................................. 62
1.5.7. Electrical contact coefficient .................................................................. 69
1.6. FORTRAN- C/C++ user subroutine interface ...................................................... 72
1.6.1. The interface ...................................................................................... 72
1.6.2. Example ............................................................................................ 75
2. User defined kinematic .............................................................................................. 82
2.1. User defined kinematic v1 ............................................................................... 83
2.1.1. Graphical interface .............................................................................. 84
2.1.2. Tags and attributes .............................................................................. 89
2.1.3. The .kin Data ...................................................................................... 92
2.1.4. Example of the .kin data ....................................................................... 94
2.1.5. Solver and subroutine ........................................................................... 98
2.1.6. Example ........................................................................................... 110
2.2. Translation files ........................................................................................... 121
2.2.1. Language dependent content ................................................................ 121
2.2.2. Creating translation files ..................................................................... 122
2.3. User defined kinematic v2 ............................................................................. 124
2.3.1. Graphical interface ............................................................................. 124
2.3.2. Tags and attributes ............................................................................. 126
2.3.3. The .kin data ..................................................................................... 130
2.3.4. Example of the .kin data ..................................................................... 131
2.3.5. Solver and subroutine: sfukin2 ............................................................. 132
2.3.6. Kin file: reading and saving information ................................................ 133
2.3.7. Example ........................................................................................... 134
3. User defined material .............................................................................................. 150
3.1. User-defined material law .............................................................................. 151
3.1.1. Material information ........................................................................... 152
3.1.2. User subroutine ................................................................................. 157
3.1.3. Example ........................................................................................... 158

ii
Scientific Tutorial
2023.2

1 User subroutines
2023.2 User subroutines FORTRAN guide (for FE-Solver)

Educational Objectives

Uses of user subroutine in Simufact Forming for better process control

Prerequisites

FORTRAN programming knowledge, general programming knowledge

1.1. FORTRAN guide (for FE-Solver)

Simufact Forming supports user subroutines written in the programming language FORTRAN.

Please note: In order to make use of and compile subroutines additional software is required which is not includ-
ed and not licensed as part of Simufact Forming. For example a compiler and a development environment is nec-
essary. For more detailed information on the software requirements, installation and configuration, please refer to
the Installation_guideline_Simufact_Forming_en.pdf. After installing all required external software you can test
your configuration using the subroutine compile_test.f provided in "C:\Program Files\simufact\forming\2023.2\sf-
Marc\tools\test". The usage of the subroutine is commented in the file itself. It can be used with any process.

In this manual some basic FORTRAN commands are reviewed. Also various user subroutine examples are discussed
in detail in the subsequent sections. For every example, the user subroutine source code and a working project file for
Simufact Forming can be accessed via the scientific branch in the Simufact demos browser. The examples will only
run successfully, in case all relevant external software is installed.

Figure 1.1. Simufact demos

In addition to this tutorial, relevant information can be found in the following places:

Templates for available subroutines in:

• "C:\Program Files\simufact\forming\2023.2\sfMarc\user"

Information on file units in the Marc Volume A (search for "File Units"):

2
2023.2 User subroutines Common block

• "C:\Program Files\simufact\forming\2023.2\doc\TechnicalReferences\sfMarc"

General information on subroutines in the Marc Volume D:

• "C:\Program Files\simufact\forming\2023.2\doc\TechnicalReferences\sfMarc"

Available common blocks:

• "C:\Program Files\simufact\forming\2023.2\sfMarc\common"

1.1.1. Common block

The common block in FORTRAN provides global storage space that can be used to share large amount of information
between different functions. This reduces the need to pass information through program arguments. Changing the
value of a common block variable in one subroutine changes the global value of the variable. Predefined solver specific
common blocks can be used in a subroutine by the command #include as:
#include "creeps.cmn"

There are many common blocks that are used by the solver to store runtime information. The creeps.cmn common
block, for example, consists of all the time related variables used in the simulation. All available common blocks can
be found in the Simufact Forming installation directory at:
C:\Program Files\simufact\forming\2023.2\sfMarc\common

Please note that required parameters, e.g. the current increment, current body-id or current element-number, is often
available as an argument of the used subroutine. If it is not available, you might be able to access it through a common
block. Examples for relevant variables are:

Table 1.1.
Common-Block Variable Description
creeps.cmn timinc Duration of current increment
concom.cmn inc Current increment number
prepro.cmn kou File unit for .out-file
machin.cmn iprtl File unit for .log-file

A typical common block is defined as:

C***********************************************************************
C
C File: creeps.cmn
C
C MSC.Marc include file
C
real*8 cptim,timinc,timinc_p,timinc_s,timincm,
1 timinc_a,timinc_b
integer icfte,icfst,
2 icfeq,icftm,icetem,mcreep,jcreep,icpa,icftmp,icfstr,
3 icfqcp,icfcpm,icrppr,icrcha,icpb,iicpmt,iicpa
real*8 time_beg_lcase,time_beg_inc,fractol,time_beg_pst
real*8 fraction_donn,timinc_ol2
c
integer num_creepsr,num_creepsi,num_creeps2r
parameter(num_creepsr= 7)
parameter(num_creepsi=17)
parameter(num_creeps2r=6)
common/marc_creeps/cptim,timinc,timinc_p,timinc_s,timincm,
1 timinc_a,timinc_b,icfte,icfst,
2 icfeq,icftm,icetem,mcreep,jcreep,icpa,icftmp,icfstr,

3
2023.2 User subroutines Common block

3 icfqcp,icfcpm,icrppr,icrcha,icpb,iicpmt,iicpa
common/marc_creeps2/
x time_beg_lcase,time_beg_inc,fractol,time_beg_pst,
x fraction_donn,timinc_ol2
c
c cptim Total time at begining of increment.
c timinc Incremental time for this step.
c icfte Local copy number of slopes of creep strain rate function
c versus temperature. Is -1 if exponent law used.
c icfst Local copy number of slopes of creep strain rate function
c versus equivalent stress. Is -1 if exponent law used.
c icfeq Local copy number of slopes of creep strain rate function
c versus equivalent strain. Is -1 if exponent law used.
c icftm Local copy number of slopes of creep strain rate function
c versus time. Is -1 if exponent law used.
c icetem Element number that needs to be checked for creep
convergence
c or, if negative, the number of elements that need to
c be checked. In the latter case the elements to check
c are stored in ielcp.
c mcreep Maximum nuber of iterations for explicit creep.
c jcreep Counter of number of iterations for explicit creep
c procedure. jcreep must be .le. mcreep
c icpa(1-6) Pointer to constants in creep strain rate expression.
c icftmp Pointer to temperature dependent creep strain rate data.
c icfstr Pointer to equivalent stress dependent creep strain rate
data.
c icfqcp Pointer to equivalent creep strain dependent creep strain
c rate data.
c icfcpm Pointer to equivalent creep strain rate dependent
c creep strain rate data.
c icrppr Permanent copy of icreep
c icrcha Control flag for creep convergence checking , if set to
c 1 then testing on absolute change in stress and creep
c strain, not relative testing. Input data.
c icpb(1-4) Pointer to storage of material id cross reference numbers.
c iicpmt creep law type ID
c =1 - power law
c =2 - solder
c =3 - steady-creep
c =4 - hyperbolic steady-creep
c iicpa Pointer to table IDs for constants in creep strain rate
c expression
c
c
c time_beg_lcase time at the beginning of the current load case
c time_beg_inc time at the beginning of the current increment
c fractol fraction of loadcase or increment time when we
c consider it to be finished
c time_beg_pst time corresponding to first increment to be
c read in from thermal post file for auto step
c
c timinc_old Time step of the previous increment
C
C***********************************************************************
!!

4
2023.2 User subroutines File operations

Information about some of the important variables defined in the common block are also explained through the com-
ment.

Users can also create their own local common blocks to store temporary run time information such as counter, variables
etc. This will help to preserve the stored informations when the function is revisited.

common /sfuser_1/ all integer variables

common /sfuser_2/ all real variables
common /sfuser_3/ all character variables
common /sfuser_4/ all pointer names

1.1.2. File operations

The FE-solver (marc.exe) uses various local auxiliary files for data storage. File unit numbers are used by the solver to
create and access these auxiliary files during the runtime. From 0 to 130 file units are used by the solver. For example,
unit 0 (variable iprtl) is used for the *.log-file of the job, whereas unit 6 (variable kou) is used for the *.out-file.
For detailed descriptions of these units, please refer to Marc Volume A: Theory and User Information: Program
Initiation - File Units. Users are advised to use file units greater than 500 to avoid conflicts with those predefined
by the solver.

To open a file, the file name must be defined. A file created during the simulation has the job name associated with
it. The various file names are stored in the common block jname.cmn. Variable kou (= unit 6) and iprtl (= unit 0) are
defined in the prepro.cmn and machin.cmn common blocks respectively. Consider the following code block:

#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif

#include "jname.cmn" ! variable ilen, jidnam

#include "machin.cmn" ! variable iprtl
#include "prepro.cmn" ! variable kou

character*80 filen,cline
parameter (iunit=57)

filen=jidnam(:ilen(1))//'.dbg'

c To write information in Format 1001 in the log and output file

write(iprtl,1001)filen(:ilen(1))
write(kou,1001)filen(:ilen(1))
.....
.....
1001 format(/13x,
*'##########################################################',
* /13x'Debug information for myself will be written in file',
* /13x,a,' in Format for Version 2023.2'
* /13x
*'##########################################################')

In the above code, the command format is used to format the output. The arguments are:

• 13x: 13 blank character space

• a: character, text is written inside quote ' '

• /: used for new line break

5
2023.2 User subroutines Conditional statement

Both FORTRAN 77 and 90 standards can be used in the subroutine. FORTRAN 77 requires the command to begin
with the 7th character and finish within the 72nd character from the left. If the command is longer then the continuation
character * (compare in above example) has to be put in the next line at 6th character from the left and then the command
can continue (from 7th character). To read a file, the following code-snippets can be used in the subroutine:

c To read a line in the file

open(iunit,file=filen,access='sequential',
* status='unknown',form='formatted')

c Read a line, jump to 110 when end of the file

100 read(iunit,'(a)',end=110) cline
goto 100
110 continue

Debugging

Writing into a file can be an easy way for debugging, although it can slow down the performance of your subroutine.
The file can be used, to print the current value of used variables or to print informative messages, e.g. in case a
condition is met. In order to create a file for debugging, you can specify a file unit and use the FORTRAN write
statement as follows:

open(unit=3001, name="_debug.txt")
write(3001, *) "Hello World"

where in this example 3001 is the file unit and * is used to specify that default formating should be applied. The file
is created in the _Run_-folder of the simulation.

1.1.3. Conditional statement

A conditional statement allows the program to make a decision. There are two forms: the "single if" statement and
the "block if" statement.

The syntax of a single if statement is:

if (logical expression) then

code-block
endif

Syntax of a block if statement:

if (logical expression) then

code-block
else if (logical expression) then
code-block
.....
.....
else
code-block
endif

Example of a block if statement:

if (i .lt. 5) then
write(*,*) "less than five"
else if(i .gt. 6) then
write(*,*) "more than five"
else if(i .eq. 5) then
write(*,*) "equal to five"

6
2023.2 User subroutines Loops

else
write(*,*) "invalid"
endif

1.1.4. Loops
Loops are used for repeated execution of a block of the code. In FORTRAN this can be achieved by a "do loop", a
"do-while loop" or a "if loop".

do loop

There are two methods to write a do loop (with and without a label):

do label variable = start, end, increment

code-block
label continue

do variable = start, end, increment

code-block
enddo

Examples: The following code prints all even numbers between 1 to 20 in decreasing order:

do 99 i = 20, 1, -2
write(*,*) 'i =', i
99 continue

do i = 20, 1, -2
write(*,*) 'i =', i
enddo

do-while loop

The while loop is used together with a do statement. The do-while loop executes the range of the do construct until
the condition inside the while statement remains true. The do-while loop can also be written in two ways:

while (logical expr) do

code-block
enddo

do while (logical expr)

code-block
enddo

Example: The code prints all even numbers between 1 to 20 in decreasing order:

i=20
do while (i.gt.1)
write(*,*) 'i =', i
i = i-2
enddo

if loop

A loop can also be written through an if-goto statement.

label if (logical expr) then

code-block

7
2023.2 User subroutines Domain decomposition method
(DDM)

goto label
endif

Example: The code prints all even numbers between 1 to 20 in decreasing order:

i=20
99 if (i.gt.1) then
write(*,*) 'i =', i
i = i-2
goto 99
endif

1.1.5. Domain decomposition method (DDM)

Simufact Forming uses the domain decomposition method (DDM) as a parallelization to reduce the simulation time.
DDM is a parallel and scalable technique in which the system first divides the problem, that can be solved separately
and at the end merges the calculated results. DDM simulations can be used to shorten the simulation time by running
the simulation on multiple processors. Simufact Forming supports the parallel execution of a simulation through DDM,
only limited by the availability of processing units in the system. The Intel MPI service should be running on the
system for using DDM for Simufact Forming simulations. The DDM can be activated in the forming control under the
parallelization tab. Each domain can then further parallelize through the shared memory parallelization (SMP) option.
The number of cores defined in this option will be allocated for each domain.

When using DDM in combination with a subroutine you have to make sure that the implications mentioned
below are considered in your subroutine.

Figure 1.2. Activation of parallel computing (DDM) in the forming control

Silent features of using DDM with the user subroutine are:

• Each domain is a separate job and runs independently to the other domains. Memory is not shared among the
domains.

• User subroutines can be used in the same way as without DDM with special care for user-defined local file handling.

• The first domain job is labeled as a parent (master) job. All domains will create runtime files with a domain number
as prefix in the job name.

For a DDM simulation, important domain related informations can be accessed in the user subroutine through the
following variables:

nprocd (from cdominfio.cmn) = number of processor used by the job

> 0 it is a DDM job
iparent (from cdominfo.cmn) = 0, then it is a parent job, important if an
activity has to be done by one domain only
iprcnm (from cdominfo.cmn) = gives the domain number in case of a DDM job
iddmarc (from gary.cmn) = flag for DDM with ARC file control

8
2023.2 User subroutines Domain decomposition method
(DDM)

= 0, no such control
= 1, first run
= 2, restart run
jidnam (from jname.cmn) = jobname of each domain
jidnam _base (from jname.cmn)= base jobname

Following is the code snippet for creating the name of the *.ARC-file by the solver. Based on the job type (DDM or
non DDM), the file names can be different.

filen=jidnam(:ilen(1))//'_'//cbody(:ilenc)//
* '_'//cinc(:ilen0)//'.ARC'
if(nprocd.gt.0)then ! in case of DDM job
ilen1=last_char(jidnam_base)
filen=jidnam_base(:ilen1)//'_'//cbody(:ilenc)//
* '_'//cinc(:ilen0)//'.ARC'
endif

To exchange runtime information among the various domains, a function domflag() can be called within the sub-
routine. When a domflag() function is called in a DDM job, the calculation in all domains will be paused at the in-
stance. Values passed in domflag() arguments from all the domains will be analyzed and user-defined mathematical
operations (like maximum, minimum, average, summation etc.) will be done. The results will then be returned to all
the domains and then the calculation will continue. The call to domflag() with dummy arguments can also used to
synchronise the domain position during calculations in case one domain might take longer time than other. In this
case the faster domains will wait for the slower domain to finish the calculations before continuing. The interface of
a domflag() function is

domflag(iflag,dflag,ichck1,ichck2,num1,num2)
c***********************************************************************
c collects iflag and dflag from domains, performs mathematical
c operations within message passing, broadcast results to domains
c
c iflag = is an integer array to be shared among domains
c dflag = is a real array to be shared among domains
c ichck1 = is an integer array of same size as array iflag. Value in
this array determines which mathematical operation to
perform on integer array (iflag)
c ichck2 = is an integer array of same size as array dflag. Value in
this array determines which mathematical operation to
perform on real array (dflag)
c
c mathematical operation flags for ichck1 and ichck2 are
c ichck = 0: minimum, return the minimum value
c ichck = 1: maximum, return the maximum value
c ichck = 2: average, return the average
c ichck = 3: sum , add and return the sum
c ichck = 4: minimum above zero
c ichck = 5: sets exit message
c ichck = 6: maximum absolute value (sign is preserved)

c num1 and num2 is the size of integer and real arrays respectivly.
c
c Results of mathematical operation, can be accessed within a domain
c from the content of array iflag and dflag.
c
c***********************************************************************
if (nprocd.gt.0)call domflag(matid,dummy,1,0,1,0) !only for ddm job
if(nprocd.gt.0)then
dbuff1(1)=test_min
dbuff1(2)=test_max

9
2023.2 User subroutines Debugging

itest1(1)=0 !get minimum value from all domain

itest1(2)=1 !get maximum value from all domain
call domflag(idum,dbuff1,0,itest1,0,2)
test_min=dbuff1(1)
test_max=dbuff1(2)
endif

The process flow path for a DDM job during the simulation is shown in the figure below:

Figure 1.3. Data flow in a DDM job

The parts of the data flow are:

• exeddm: A master program that controls the overall data flow

• marc: FEM solver runs analysis in DDM fashion and returns to sf_exeddm when remeshing is required

• sfremap: Program that merges domain-based *ARC-files, calls sfmultimesh and does the data mapping, returns to
sf_exeddm with a merged *.ARC-file

1.1.6. Debugging
It is possible to locate the cause of run-time and logical errors in a user subroutine through a debugger. With the help of
a debugger, variable values and information flow through the subroutine can be monitored during the simulation. This
helps to trace the source of a logical error in the subroutine. A detailed process to set up the debugging environment
with the debugger in Microsoft Visual Studio (2022) is described in this section. Please note that only non DDM jobs
can be debugged with this method.

1. Setup the debug environment: It is recommended to create two batch files and put them in a user path directory.
This will enable to call batch files from any directory through the command line. The first is a command to invoke
the Intel FORTRAN compiler and the second is to start the debugger.

a. To setup a Intel FORTRAN compiler environment at a command line:

Create a first batch script file and name it as intel_comp.bat, add the following two lines to it and save it in
the user path directory:

set INTELPATH=C:\Program Files (x86)\Intel\oneAPI

if exist "%INTELPATH%compilervars.bat" call "%INTELPATH%compilervars.bat"
intel64 vs2022

From a command line console (cmd.exe) executing intel_comp.bat should make the Intel FORTRAN compiler
accessible at the command prompt. If not, please check for the correct installation directory of the compiler and
modify the above code. Also check if the location of intel_comp.bat file is included in the local user path profile.

b. In the second batch script, set the Simufact Forming installation path available for the debugger:

@echo off
setlocal

10
2023.2 User subroutines Debugging

call intel_comp
set FORMINGDIR=C:\Program Files\simufact\forming\2023.2
set SDIR=D:\log
set BINDIR=%FORMINGDIR%\sfMarc\bin\win64i8
set EXITMSG=%SDIR%\MESSAGES
set LIBDIR=%FORMINGDIR%\sfMarc\lib\win64i8
set CUDADIR=%FORMINGDIR%\sfMarc\lib\win64i8\cuda
set MESHERDIR=%FORMINGDIR%\sfTools\sfMeshing\bin
set PATH=%MESHERDIR%;%BINDIR%;%LIBDIR%;%SDIR%;%PATH%;%CUDADIR%
call devenv
endlocal

The above code sets up the compiler environment and adds various Simufact Forming libraries to the existing
system path. At the end it calls the Visual Studio debugger with the command devenv. Running the above script
should open the Visual Studio debugger. Modify the various paths depending upon the local installation and
settings. Give the script a suitable name and store it in a folder included in the local user path profile.

The above setup have to be done once for a system.

2. Enable the debug flag in the Simufact Forming script: By default the debug flag is switched off in Simufact
Forming script and should be manually enabled. Go to the Simufact Forming installation directory at C:\Program
Files\simufact\forming\2023.2\sfMarc\tools, open the file include_win64i8.bat and search for
the option MARCDEBUG. Remove the REM in front of the line set MARCDEBUG=ON. The modified script should
look as:

The above modification will create the solver binary with the debug flag for the compiled subroutine code. Please
note that the job will run slower when the subroutine is compiled with the debug flag. It is recommended to switch
off the debug flag for actual simulations.

3. Create the debug executable:

a. Setup the simulation job in Simufact Forming. Attach the user subroutine that needs to be debugged in the solver
option in the forming control.

b. Create the complete input files in the directory _Run_ of the project by submitting the job. Various job files are
created in the _Run_ directory. The simulation should then be safely stopped.

c. The _Run_ directory should contain the job input file (*.dat), the material file (*.umt), a copy of the user
subroutine and two batch files named runsf.bat and cleansf.bat. Open the file runsf.bat in a text editor and
copy and modify the script as:

call "C:\Program Files\simufact\forming\2023.2\sfMarc\tools\run_marc.bat"

-u "ufric" -j friction -b no -save yes

The last two arguments -b no -save yes are added arguments. The argument -save yes prevents the
deletion of the executable by the Simufact Forming script after the job is completed.

d. Open a command console (cmd.exe) in the _Run_ directory and execute the above modified code in the com-
mand prompt. This should create a solver *.exe together with debug files with the extension *.pdb. The simula-
tion can be terminated once these files are created. Double-clicking on cleansf.bat will clean the _Run_ directory
with the temporary job runtime log/output files.

4. Debugging in Visual Studio debugger:

a. Open the Visual Studio debugger through the script described in the first step.

b. Open the created solver executable (*.exe) in Visual Studio as Project/Solution as:

11
2023.2 User subroutines Variable units

c. Open the subroutine (*.f) as file:

d. Add the job file name (*.dat) as project argument with the argument -j in the project properties as

5. Now the project is ready for debugging in Visual Studio. Insert break points at different sections of the user sub-
routine that need to be analysed. The debugger will pause when that portion of the subroutine is executed.

It may happen that Visual Studio complains of some missing *.dll-files or libraries. Check that the necessary runtime
dependencies/distributables are installed. If the error complains of missing functions in existing libraries then it is due
to the libraries in the system are old. Updating the program with the latest version that initially provided the libraries
should get rid of the error.

1.2. Variable units

By default Simufact Forming uses SI units for variables in the simulation. The use of incorrect units in the simulation
can result in wrong simulation results. This becomes important when input values are manually changed in the input
file of the job (*.dat) or in the material file (*.umt) before the simulation is started or user subroutines are used
to perform the calculation. Please note that Simufact Forming and SimufactGP uses different default variable units.
Hence, if a subroutine is shared between Simufact Forming and SimufactGP, special care should be taken for the
variable units. The default units of the quantities used in the simulation are tabulated below.

Variable units

12
2023.2 User subroutines Coordinate axis

Quantity SI unit mm/tonne/s/K or SI-mm Conversion factor

(Simufact Forming) (SimufactGP) from SI to SI-mm
Length m mm x 103
Time s s x1
Mass kg tonne = Mg x 10-3
Force kg·m/s2 = N tonne·mm/s2 = N x1
3 3
Density kg/m tonne/mm x 10-12
Stress kg/m/s2 = N/m2 tonne/mm/s2 = N/mm2 x 10-6
Energy kg·m2/s2 = J tonne·mm2/s2 = mJ x 103
Temperature K °C - 273.15
2 2 2 2
Specific heat capacity m /s /K mm /s /°C x 106
Heat convection kg/s3/K tonne/s3/°C x 10-3
Thermal conductivity kg·m/s3/K tonne*mm/s3/°C x1
Thermal expansion m/m/K mm/mm/°C x1

1.3. Coordinate axis

2D simulation in Simufact graphical user interface (GUI) is set up in z-x plane. However, solver does the calculation
in x-y plane. The output generated by the solver is transferred back to z-x plane in Simufact GUI. In case user extract
post variables (for example coordinates) for a 2D simulation through user subroutine, the solver will return the x and
y component but not the z component. The returned x and y components are in fact z and x component with reference
to Simufact GUI.

Job type Forming GUI Solver

2D

3D

No such transformation occurs in case of 3D jobs.

1.4. User subroutine and utility routines

User subroutines allow flexibility to customize the simulation process for non-standard problems. Users can either
introduce new physics in the simulation or can get customized process output according to the need. User subroutines

13
2023.2 User subroutines Subroutine details

are to be written with FORTRAN 77 or 90 standard, that can be easily linked with the solver to bypass the inbuilt routine
with the proper control setup. A comprehensive list of all the available subroutine files and its structure can be found
in the Marc manual "Volume D", see "C:\Program Files\simufact\forming\2023.2\doc\TechnicalReferences\sfMarc" .
Interfaces of these subroutines are compatible when used in a Simufact Forming simulation. Subroutines written in
C/C++ are also supported if the functions are embedded inside a DLL (Dynamic Link Library). For more details refer
to the section "FORTRAN-C/C++ user subroutine Interface" of this tutorial.

Special utility routines are available to retrieve nodal and elemental process variables during the simulation. These
routines are extremely helpful to retrieve process information in user subroutines as they may not be available through
subroutine's interface.

1.4.1. Subroutine details

To use a subroutine successfully, it is important that the correct version of Intel FORTRAN Compiler is installed on
the system. On Windows, Visual Studio is additionally required to provide linker. Linker attach the compiled object
file of the subroutine to the standard solver. Consider the original user subroutine ueloop.f from Marc, which is called
over all the elements in a loop. This subroutine can be used to define or modify variables data stored in the common
blocks during various elemental loop.

subroutine ueloop(m,n,iflag)
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
c ** Start of generated type statements **
integer iflag, m, n
c ** End of generated type statements **
c
c user routine that gets called in key element loops
c may be used in a variety of ways
c
c m element number
c n internal element number
c iflag 1 - opress loop
c 2 - oasemb loop
c 3 - oasmas loop
c 4 - ogetst loop
c 5 - scimp loop
c
return
end

The subroutine is divided into different sections. The top part contains the call arguments and the common block. Call
arguments describe the list of variables that are passed to the subroutine, for example m, n, iflag are variables that
are passed as arguments in this subroutine:

subroutine ueloop(m,n,iflag)

Common blocks are used to access the required variables as it is not possible to pass all the variables through a call
argument. Almost all information is available through the usage of proper common blocks. Common blocks can be
activated by adding an #include command in the user subroutine. The syntax to use in the user subroutine is:
#include "yyy.cmn", where yyy.cmn is the name of the common block as:

#include "creeps.cmn"

A list of all common blocks used internally by the solver can be found at C:\Program Files\simufact\form-
ing\2023.2\sfMarc\common. For more details about common blocks, refer to the Marc manual "Volume D".

14
2023.2 User subroutines Subroutine details

Comments in the subroutine help to explain the variables and program functionality for easier understanding and
readability of the program code. Comment lines are beginning with a character c. Users can add their own code before
the return statement. An example of the extended user subroutine to calculate total elastic strain is given below. Two
subroutines are used to get the desired results. The first subroutine (ueloop) calculates and stores the incremental
plastic strain (eplas) and the other (plotv) outputs the values to the result for visualization.

subroutine ueloop(m1,m2,iflag)
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
c ** Start of generated type statements **
#include "hards.cmn" !nstats
#include "creeps.cmn" !timinc
#include "arrays.cmn" !iintel
#include "elemdata.cmn" !ieltype
#include "space.cmn" !ints
#include "concom.cmn" !ipass
c
integer iflag, m1, m2
c ** End of generated type statements **
integer layer, icode,ioflag,ibodyid,nn1,inn,i1,ityp,
* istat
real*8 eplas,deplas
c
c user routine that gets called in key element loops
c may be used in a variety of ways
c
c m1 user element number
c m2 internal element number
c iflag 1 - opress loop
c 2 - oasemb loop
c 3 - oasmas loop
c 4 - ogetst loop
c 5 - scimp loop (end of increment)
c
c IPASS = 1 STRESS
c IPASS = 2 HEAT TRANSFER
c IPASS = 3 FLUIDS
c IPASS = 4 JOULE HEATING
c IPASS = 5 DIFFUSION
c IPASS = 6 ELECTROSTATICS
c IPASS = 7 MAGNETOSTATIC
c IPASS = 8 ELECTROMAGNETICS

if(ipass.ne.1) goto 999 !only for stress pass

if(iflag.eq.5.and.nstats.gt.1)then
c utility routine to extract contact body ID
c
c ielem element number
c ibodyid body ID
c ioflag =0, user element number
c =1, element number is internal
ioflag=0
call get_bodyid(m1,ibodyid,ioflag)
if(ibodyid.ne.1)goto 999
c

15
2023.2 User subroutines Subroutine details

c get number of integrationpoints for element (nn1)

ityp=ieltype(m2)
i1=ityp-1
nn1=ints(i1+iintel)
c
do inn=1,nn1
layer=1
icode=28 !strain rate
call get_elmvar(icode,m1,m2,inn,layer,deplas)
c calculate incremental plastic strain
deplas=deplas*timinc
c
c read data ioflag = 0
c write data ioflag = 1
c
c get 2nd state variable (ioflag=0)
istat=2
icode=1000000+istat
ioflag=0
call elmvar1(icode,m1,inn,layer,eplas,ioflag)
c cummulate value
eplas=eplas+deplas
c store 2nd state variable (ioflag=1)
istat=2
icode=1000000+istat
ioflag=1
call elmvar1(icode,m1,inn,layer,eplas,ioflag)
enddo
endif

999 return
end

subroutine plotv(v,s,sp,etot,eplas,ecreep,t,m,nn,kcus,ndi,
* nshear,jpltcd)
c* * * * * *
c
c define a variable for contour plotting (user subroutine).
c
c v variable to be put onto the post file
c s (idss) stress array
c sp stresses in preferred direction
c etot total strain (generalized)
c eplas total plastic strain
c ecreep total creep strain
c t array of state variable (temperature first)
c m(1) user element number
c m(2) internal element number
c m(3) material id
c m(4) internal material id
c nn integration point number
c kcus(1) layer number
c kcus(2) internal layer number
c ndi number of direct stress components
c nshear number of shear stress components
c jpltcd the absolute value of the user's entered post code
c

16
2023.2 User subroutines Subroutine details

c
c the components of s, sp, etot, eplas and ecreep are given in the
order
c as defined in volume B for the current element.
c
c* * * * * *
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif

c ** Start of generated type statements **

real*8 ecreep, eplas, etot
integer jpltcd, kcus, m, ndi, nn, nshear
real*8 s, sp, t, v
c ** End of generated type statements **
dimension s(*),etot(*),eplas(*),ecreep(*),sp(*),m(4),kcus(2),
* t(*)
integer ioflag,icode

if(jpltcd.eq.1)then
c get elastic strain
icode=127 !equivalent elastic strain post code
ioflag=0 !read
call elmvar1(icode,m(1),nn,kcus(1),v,ioflag)
v=v+t(2) !add 2nd saved state variable
endif
return
end

Several Simufact Forming specific subroutines have been developed. sfukin.f, user_fstress.f and urigid_rotate.f are
few among them. sfukin.f can be used for user kinematics as discussed in Chapter 2.

The subroutine user_fstress.f can be used to compute the flow stress based on a user equation. This subroutine supports
the material input parameters from Simufact Material and Simufact Forming.

subroutine user_fstress(ilst,ityp1,t1,ep1,er1,p1,fstress1)
c
c subroutine to compute flow stress based on user equation
c ilst(1) - material id
c ilst(2) - element id
c ilst(3) - integration point
c ilst(4) - layer id
c ilst(5) - internal layer number
c ityp1 - equation type (>100)
c Note
c This is the "Model number" defined in the GUI + 100
c E.g. If "2" is defined in the GUI, ityp1 is 102
c t1 - temperature
c ep1 - plastic strain
c er1 - plastic strain rate
c p1 - material parameters from input
c fstress1 - flow stress

The subroutine urigid_rotate.f can be used to define rotational motion during the simulation.

urigid_rotate(phix,phiy,phiz,raxis,rangle,
& rsystem,timinc,ctime,idie,node,iflag)

17
2023.2 User subroutines Subroutine integration

c Input:
c idie = load controlled die id
c node = load controlled die - rotational control node id
c iflag = 0: compute rotation angle and axis
c = 1: compute only the local system rotation matrix

c Output:
c raxis(3) = rotation axis
c rsystem(3,3) = local system rotation matrix

1.4.2. Subroutine integration

This section explains the method to integrate the subroutine in the simulation process. For this, the subroutine has
to be written first. All available dummy subroutines are located in the installation folder of Simufact Forming at C:
\Program Files\simufact\forming\2023.2\sfMarc\user. It is recommended that before modifying
any subroutine, it should be copied to another working directory, so that the original copy of the file still exists. The
copied subroutine file can be linked with the process in the forming control under Advanced > User-defined:

Figure 1.4. Linking a user subroutine in the forming control

If the FORTRAN user subroutine is copied to the directory AdditionalFiles of the project, then just writing the complete
file name of the subroutine in the forming control is sufficient. The Simufact Forming graphical user interface (GUI)
will search in the directory AdditionalFiles of the project first.

Figure 1.5. Subroutine in directory AdditionalFiles

18
2023.2 User subroutines User variables, state variables and post
variables

With this convention, the subroutine always remains with the project in case of transferring it to another computer.
If the complete path or file name of the linked subroutine is wrong then the text field will be highlighted in orange.
Alternatively, a pre-compiled solver (with extension .exe) can also be used. This can be linked in the forming control
window, too. A compiled solver is created in the run directory of the job when a job is run using subroutines or
alternatively can be created from the command line. This solver contains instructions from the subroutine and can be
used on machines without a FORTRAN compiler.

Figure 1.6. Linking a pre-compiled solver in the forming control

In addition to linking the subroutine file as above, certain subroutines require additional steps of activation from the
menu Advanced > User-defined in the forming control. When using the following subroutines: umotion, uhtcoe,
uhtcon, ufric, sepstr and sepfor, additionally have to be activated in the forming control under Advanced > User-
defined > User subroutine.

Figure 1.7. Activating a selected subroutine in the forming control

1.4.3. User variables, state variables and post variables

In the forming control window under Advanced > User-defined, there are options to define different types of variables
for the simulation. These variables can then be accessed through the user subroutines during runtime. There are three
types of variables: user variables, state variables and post variables. These types are explained in detail below,
starting with a short overview.

User variables are used to provide numerical input for the simulation through the GUI. This avoids manual adjustment
of the source code, when the value of the input variables changes. A possible example could be a diameter which
changes from one simulation to another. In that case a variable for this diameter can be defined in the user variables
section in the GUI. The user can change this value before starting the simulation. In the subroutine this value is
accessible through the common block "sf_userdata.cmn".

19
2023.2 User subroutines User variables, state variables and post
variables

State variables, in contrast to user variables which allow external change of one single value, can be changed only
directly through the subroutine but not externally from the GUI. However, their usage and number has to be declared in
the GUI, to make sure enough memory is allocated. State variables are used to store intermediary results or calculated
variables for each element or node at a given moment and therefore capture the state of the body at that moment. A
typical application could be a constant value which is added to each node or a user specific formula which is evaluated
for each node/element. The concept of state variables is a general concept which is also used in the solver itself, to
store the internal marc specific results. The values of these internal state variables are accessible from the subroutine
and can be used for the calculation of own result values, which can then be stored in own state variables. Please refer
to Section Post codes for more details on how to access the solver internal result variables.

Post variables can be used to display any arithmetical output in the GUI. In a typical use case the value of the state
variable is not only relevant for calculation inside the subroutine but should also be accessible as a user specific result
value in the GUI. In that case, the values of the original state variable are copied into a post variable solely available
within the subroutines upstno.f (nodal) and plotv.f (elemental). The typical work flow is as follows: The results
created in the main subroutine, e.g. ufninc (nodal) or ueloop (elemental) are written to the state variable using the
utility function ndstatvar() (nodal) or elmvar1() (elemental). The value of these state variables are then accessed from
inside the subroutine upstno (nodal) or plotv.f (elemental) using the same utility functions ndstatvar() or elmvar1() but
this time with reading flag. The obtained results are then assigned to the post variable for the solver to write into the
result file which gets displayed by the GUI. The subroutines upstno.f and plotv.f are called as often as the total number
of post variables requires. For example, for a total of 2 node variables, during the first call the result are written in
the first nodal variable and during the second call the result are written into the second nodal variable. To define the
number of required calls, the use of post variables has to be declared in the post variable section in the GUI. Also
meaningful names for each variable can be given in the post variable section. For more details, please see the example
in the section User subroutine - examples.

Variable definition in the forming control

As shown in the figure below, all three types of variables: user variables, state variables and post variables can be
defined in the Simufact GUI. Details on the programmatic use of these variables can be found below. The numbers
in the list refers to the red numbers in the figure:

Figure 1.8. Variable definition in the forming control

1. User variables: Integer, real and string variables can be defined here. The solver allocates memory based on the
number of variables defined here.

• To access these variables in a subroutine, include the common block: #include "sf_userdata.cmn"

20
2023.2 User subroutines Utility routines

• Through the above common block one can determine in the subroutine if user variables are initialized or not. If
initialized then it can be defined how many of each kind are to be initialized through:

• niuser: number of integer user variables

• nruser: number of real user variables

• nsuser: number of string user variables

• The content of the variables can then be accessed in subroutines as:

• iuserdata(i): provide content of ith integer data

• ruserdata(i): provide content of ith real data

• suserdata(i): provide content of ith string data. At present the length of the character array is set to 80.

2. State variables: The solver allocates extra memory to be used for information storage during the simulation based
on the number defined here. The state variable values can be read and written through the utility routines and post
codes. At a remeshing step the values stored in the state variables will also be mapped with the same mapping
function as other simulation variables like stress, strain etc. The values stored in these variables are of real type.
State variables can be defined for nodes as well as for elements:

• Nodal state variables for extra memory at each nodes. The variable nstatev from sf_userdata.cmn stores the
number of user-defined nodal state variables.

• Elemental state variables for extra memory at each integration point of the elements. By default one state variable
is always initialized and is used to store the temperature. So, if the user needs 2 state variables, the value should
be set to 3. The variable nstats from hards.cmn stores the number of user-defined elemental state variables.

Relevant utility routines for accessing and writing into user-defined state variables are ndstatvar (nodal) and elm-
var1 (elemental). Please refer to the section Utility routines for more details on access of state variables and ex-
amples.

3. Post variables: These variables trigger the output of user-defined post variables through the subroutines upstno.f
(nodal) and/or plotv.f (elemental) output. At present only 20 post variables of either kind is supported from the
Simufact Forming interface. Refer to section Post codes for common post code definition.

1.4.4. Utility routines

Following is the list of Simufact Forming specific utility routines through which all the nodal and elemental process
variables can be extracted. The examples provided in the next chapter uses these utility routines.

1.4.4.1. elmvar1()
elmvar1 should be used to extract elemental values from the simulation. It can also be used to read and write the user-
defined elemental variables. This can also be used to extract principle stress and strain at the element level.

elmvar1(icode,mmx,nny,kkcy1,var,ioflag)

Input:
icode = post code
mmx = user element number
nny = integration point number
kkcy1 = layer number
ioflag = 0, read data
= 1, write data

21
2023.2 User subroutines Utility routines

Output:
var = output or save value, depending upon ioflag

userdefined elemental post code starts from 1000000 for example:

icode = 1000000+1 => temperature
= 1000000+2 => first user defined elemental variable and so on

1.4.4.2. get_elmvar()
get_elmvar should be used to get element values at the integration points.

get_elmvar(icode,mmx,mmy,nny,kkcy1,var)

Input:
icode = post code
= 0 (pick up integration point coordinates)
mmx = user element number
mmy = internal element number
nny = integration point number
kkcy1 = layer number

Output:
var = output value(s)

1.4.4.3. ndstatvar()
ndstatvar should be used to read and write the user-defined nodal variables.

ndstatvar(istat,ind,var,ioflag)

Input:
istat = user defined nodal variable id
= 1 ,first user defined node variable
ind = node number
ioflag = 0, read data
= 1, write data

Output:
var = output or save value, depending upon ioflag

1.4.4.4. get_nodvar()
get_nodvar should be used to extract nodal values based on the user node numbers, nodal post codes and/or element
post codes. For element post codes, the integration point values will be extrapolated to the nodes and averaged.

get_nodvar(lnode,lnpost,lepost,var_node,var_elem,num_node,
num_npost,num_epost,nsize_node,iflag)

Input:
lnode(num_node) = array of node numbers for which nodal and elemental
values are to be extracted
lnpost(num_npost) = array of nodal post codes
lepost(num_epost) = array of elemental post codes
num_node = number of nodes for which post code is requested
num_npost = number of nodal post codes requested for each node
num_epost = number of element post codes requested for each node
nsize_node = total dimension of nodal post results

22
2023.2 User subroutines Utility routines

iflag = 0 node numbers are external

= 1 node numbers are internal

Output:
var_node(nsize_node,num_node) = return nodal results from nodal post codes
var_elem(num_epost,num_node) = return nodal results from element post
codes

1.4.4.5. get_bodyid()
This utility routine extracts the contact body ID. It is important to know if the element belongs to the workpiece or to
other bodies. In a single workpiece job, the workpiece elements will have the body id 1.

get_bodyid(ielem,ibodyid,iflag)

Input:
ielem = element number
iflag = 0, element number is external
= 1, element number is internal

Output:
ibodyid = body ID

1.4.4.6. get_material()
This utility routine extracts the material ID and number.

get_material(ielem,matnum,matid,iflag)

Input:
ielem = element number
iflag = 0, element number is external
= 1, element number is internal

Output:
matnum = material number
matid = material id

1.4.4.7. read_kin_version2()
read_kin_version2 is used to read the user kinematic kin file format, version 2. This routine will read rest of the kin
file and copy the information in modules to be available globally. The routine only reads second and third data block
entry. The first block is read elsewhere, so it is important to call this routine only in sfukin2 subroutine when called
with iflag=0.

read_kin_version2(ifunit,ierrcode)

Input:
ifunit = unit number of already opened kin-file

Output:
ierrcode = non zero if solver gets error reading kin file

1.4.4.8. index_info()
This integer function returns an integer value associated with the unique name used in kin file format 2. The search
is case in-sensitive.

23
2023.2 User subroutines Utility routines

index_info(uname)

Input:
uname = unique identifier word

Output:
index_info = index value where "uname" information is stored
returns -1 when no "uname" is found

The function returns -1 if the search fail.

1.4.4.9. get_body_bc()
This utility routine is used to get point load and fixed displacement boundary condition information for a body. Through
this information, a boundary condition of a body can be enabled and disabled on demand. Syntax of a routine is:

get_body_bc(ibd,ibc_ptld,ibc_dit)

Input:
ibd = body number

Output:
ibc_ptld = point load boundary condition index number
ibc_dit = fixed displacement boundary condition index number

The extracted boundary condition index can be saved for enabling it later, for example in force controlled bodies.

1.4.4.10. change_lcase_title()
change_lcase_title can be used to change the name of load case that appears in the .sts file.

change_lcase_title(lctitle)

Input:
lctitle = New load case title, maximum of 24 character length

1.4.4.11. sfkinbd()
To get bounding box, minimum and maximum coordinates of a body in x, y, and z direction, this utility routine can
be used. The syntax is:

sfkinbd(ibody1,xmin,xmax,iflag)

Input:
ibody1 = body id whose current min-max need to be calculated
iflag = not used. put 0

Output:
xmin(3) = minimum coordinate of body in x, y and z direction
xmax(3) = maximum coordinate of body in x, y and z direction

1.4.4.12. apply_fixed_bc()
Utility routine to apply fixed boundary condition. The syntax is:

apply_fixed_bc(ibd,ipos,vel)

24
2023.2 User subroutines Utility routines

Input:
ibd = body number
ipos = position in xkinbcvel
vel(3) = velocity (in x, y ,z)

Above ipos refer to position defined in kinematic boundary condition at start.

1.4.4.13. activate_force_bc()
This logical function is used to make a body force controlled when force switch criteria are satisfied. This routine is
a logical function which returns true when the boundary condition is successfully activated. For a body to be force
activated, point load boundary condition and fixed displacement boundary condition id is required.

activate_force_bc(ibc_ptld,ibc_dit,kdir)

Input:
ibc_ptld = ptld boundary condition index for the body
ibc_dit = dit boundary condition index for the body
kdir = direction of activation (1=x-axis; 2=y-axis, 3=z-axis)

Output:
activate_force_bc = true if activation succeeds else false

1.4.4.14. calc_ring_ypos()
calc_ring_ypos can be used to determine ring thickness along y-axis (x=0) at certain height. This logical routine
search the intercept of ring mesh along:

• y axis along x=0 plane

• at prescribed height 'h' in z direction

Figure 1.9. Ring intercept

25
2023.2 User subroutines Utility routines

Syntax of the routine is

calc_ring_ypos(atheight,yposition,ibodyid)

Input:
ibodyid = bodyid of ring
atheight = relative height in z direction where the y position will be
measured

Output:
yposition (4) = is a real array of size 4
returns the ring absolute y coordinate from left
calc_ring_ypos = TRUE if all the four points could be found else FALSE

The accuracy of result depend upon the mesh size of the ring. With this routine, ring diameter (inner and outer), left
and right ring thickness can be calculated.

1.4.4.15. new_position_on_arc()
This routine calculates the position of body when moved on an arc path with an angle.

new_position_on_arc(x0,y0,rad,x1,y1,beta,x2,y2)

Input:
x0,y0 = arc center coordinate (x,y)
rad = arc radius
x1,y1 = current position on arc
beta = move the point by ange beta (radian)

Output:
x2,y2 = new position(x,y)

1.4.4.16. twocirc_intersec()
This routine returns the intercept coordinates of two circle. The syntax of the routine is

twocirc_intersec(x1,y1,r1,x2,y2,r2,sx1,sy1,sx2,sy2,nintersec)

Input:
x1,y1,r1 = center x-, y-coordinate and radius of circle 1
x2,y2,r2 = center x-, y-coordinate and radius of circle 2

Output:
sx1,sy1 = x- and y-coordinates of intersection 1 (if existent)
sx2,sy2 = x- and y-coordinates of intersection 2 (if existent)
nintersec = number of intersections

1.4.4.17. get_val_table()
get_val_table can be used to evaluate table information that is entered through kinematic xml dialog. GUI writes
the table information in the dat file and the table id in the kin file. With table id and the value of first column, the
corresponding 2nd column value can be evaluated.

get_val_table(id_tab,col1,col2)

Input:
id_tab = table id

26
2023.2 User subroutines Utility routines

col1 = first column value

Output:
col2 = corresponding second column value

Solver will linearly interpolate the value between two points. Note: Only one independent variable for table is sup-
ported at present.

1.4.4.18. tabva2()
This routine can help to extract values from the table driven input inside the user subroutine. This routine is useful if
the user wants to reference a value from a table that is being used by the solver.

call tabva2(refval, evalue, idtable, 0, 0)

Input:
refval = reference value by which the output will be scaled
idtable = table id in the input file

Output:
evalue = evaluated value

The solver already has the temperature information in the variable tempi. If the user wants to evaluate the value at
a different temperature then the user has to define a different temperature as following:

1. Include the ctable common block #include "ctable.cmn"

2. Set the temperature in the variable tempi as tempi=user_temperature

3. Call tabva2().

If the value of tempi is changed in the subroutine, then it should be reverted back to its original value.

1.4.4.19. trackmypoint()
trackmypoint can be used to extract particle information

trackmypoint(ipointid,ipostcode,results,iflag)

Input:
ipointid = particle id
ipostcode = element postcode of result variable

Output:
results = return results() array
iflag = returns extraction status as:
Positive: Successful, size of resultarray
Negative: extraction unsuccessful
-1: Something unknown went wrong
-2: Invalid point id, point id does not exist
-3: Post code not supported
-4: Point lies inside inactive element, nothing will be
calculated

Example:

For position call the routine with ipostcode =0. Coordinates (x,y,z) will be returned in xcoords(1:3) array for point
i. Value of iflag will be set to 3.

27
2023.2 User subroutines Common user subroutines

call trackmypoint(i,0,xcoords,iflag)

This routine calculates one post code at a time.

Please note that this routine is currently not supported for DDM.

1.4.5. Common user subroutines

Some of the commonly used subroutines are described in this section.

1.4.5.1. Incremental subroutines

These subroutines are called once every increment at different stages:

• ubginc() : called at the beginning of every increment.

• ubgpass() : also called at the beginning of every increment but after calling ubginc(). This routine is called for every
pass of analysis and multiple times for a coupled analysis. A regular forming job has coupled thermal and mechanical
passes where the first call is a thermal pass with flag ipass=2 followed by a mechanical with ipass=1. So care
should be taken for which pass the calculation should be done.

• ufninc() : called at the end of every increment but before the results are created for visualization.

• uedinc() : also called at the end of every increment but after the results are created for visualization.

1.4.5.2. Elemental and nodal subroutines

These subroutines are called within an increment for every node or element accounted in the simulation.

• ueloop() : called for every element in the simulation. Also called at different stages based on the value of iflag.
With iflag=1 the routine is called after ubginc(). And before ubgpass() with iflag=5, it is called before the
results are created.

• plotv() : called for every integration point. Used to generate custom result output.

• upstno() : called for every node. Used to generate custom nodal result output.

Users can also create local elemental/nodal loops within any subroutine. Care should be taken to filter out the elements
and/or nodes of the workpiece from those of the dies. Also refer to Marc Volume D: User Subroutines and Special
Routines - Output Quantities User Subroutines for more details about postprocessing. The flowchart below shows the
schematic of the various routine calls. This is not the complete flowchart of the simulation.

28
2023.2 User subroutines Common user subroutines

Figure 1.10. Flowchart highlighting subroutine call

29
2023.2 User subroutines Post codes

1.4.6. Post codes

Post codes are used in the utility routines to get the required nodal or elemental quantities from the system. A post code
is an integer value which uniquely identifies each available variable. A comprehensive list of available result values
and the corresponding post codes can be found in Marc Volume C in "C:\Program Files\simufact\forming\2023.2\doc
\TechnicalReferences\sfMarc" or via Simufact Forming GUI -> Help -> Technical references. Examples of relevant
post codes are also listed at the end of this section.

There exist elemental as well as nodal process variables.They can be accessed and used for calculation of own result
values, see section User variables, state variables and post variables. The access to these internal state variables is
read-only. For user-defined state variables, read and write access is available. A short overview on how to access state
variables is provided below. Detailed explanations of the functions can be found in the section Utility routines.

Table 1.2. Utility routines for accessing elemental post values

Function Usage Post codes Type of Input Return Value
per call post value
get_elmvar Read internal only single post only element both internal and scalar or vector
Marc Post Codes value at once external element
number
elmvar1 Read internal only single post only element external element scalar or vector
Marc Post Codes / value at once number (internal) / scalar
Read and write (user)
user defined
state variables

Table 1.3. Utility routines for accessing nodal post values

Function Usage Post codes Type of Input Return value
per call post value
nodvar Read internal only single post only nodal external node ids scalar or vector
Marc Post Codes value at once
get_nodvar Read internal multiple post val- nodal and ele- external or inter- scalar or vector
Marc Post Codes ues at once ment (extrapolat- nal node id (based
ed to node) on value of iflag)
ndstatvar Read and write only single post only nodal internal node ids scalar
user defined value at once
state variables

Internal and external numberings

There exist two different numberings for node and element variables: Internal (marc) and external (user) node/element
numbering. The external node/element ids are based on the actual user input (GUI) and hence can be found in the .dat-
file. The internal node/element ids are only available for the solver. They are created in the process of optimizing
the formulation of the underlying mathematical problem for numerical reasons. The external ids are unique, which is
especially relevant for DDM Jobs. The numbering for internal Ids is based on the domains, i.e. the same internal id
will exist in different domains and is therefore not unique.

Depending on the purpose and application of available utility functions, either internal or external are required as
input. During programming, it is important to make sure, which type of node/element id is available (e.g. passed
as an argument to the user defined subroutine) and which type is required by the utility routines used for querying
result values. Using the wrong type of id may result in wrong calculation results and solver crashes. For converting
internal into external node ids and vice versa the integer functions nodext() and nodint() are available. For converting
element ids, ielext() and ielint() can be used, see examples below. For more details refer to Marc Volume D: Chapter
1 Introduction in "C:\Program Files\simufact\forming\2023.2\doc\TechnicalReferences\sfMarc"

30
2023.2 User subroutines Post codes

c Example for converting an external node id to an internal node id

c Conversion has to be executed for each node (loop)
integer nodint
internalNodeId = nodint(externalNodeId)

c Example for converting an internal element id to an external element number

c Conversion has to be executed for each element (loop)
integer ielext
externalElemNum = ielext(internalElemNum)

1.4.6.1. Element post codes

Some of the common used element post codes are tabulated below. For a list of all available post codes, please refer
to Marc Volume C: Model Definition Options - Program Control - POST (Model Definition) in "C:\Program Files
\simufact\forming\2023.2\doc\TechnicalReferences\sfMarc".

Table 1.4. Selected list of available element post codes

icode Description
1-6 Components of strain. For rigid-perfectly plastic flow problems, compo-
nents of strain rate
7 Equivalent plastic strain (integral of equivalent plastic strain rate). For rigid-
perfectly plastic flow problems, equivalent plastic strain rate
9 Total temperature
10 Increment of temperature
11-16 Components of stress
17 Equivalent von Mises stress
18 Mean normal stress (tensile positive)
19 User-defined variable via the plotv user subroutine. See Marc Volume D:
User Subroutines and Special Routines.
21-26 Components of plastic strain
27 Equivalent plastic strain
28 Plastic strain rate
31-36 Components of creep strain
37 Equivalent creep strain
51-56 Real components of harmonic stress
57 Equivalent real harmonic stress
121-126 Elastic strain
127 Equivalent elastic strain

1.4.6.2. Node post codes

Following are the node post codes available for nodal values. For a list of all available post codes, please refer to
Marc Volume C: Model Definition Options - Program Control - POST (Model Definition) in "C:\Program Files\sim-
ufact\forming\2023.2\doc\TechnicalReferences\sfMarc".

Table 1.5. Selected list of available node post codes

0 = Coordinates 1 = Displacement 2 = Rotation
3 = External force 4 = External moment 5 = Reaction force

31
2023.2 User subroutines User subroutine: Examples

6 = Reaction moment 7 = Fluid velocity 8 = Fluid pressure

9 = External fluid force 10 = Reaction fluid force 11 = Sound pressure
12 = External sound source 13 = Reaction sound source 14 = Temperature
15 = External heat flux 16 = Reaction heat flux 17 = Electric potential
18 = External electric charge 19 = Reaction electric charge 20 = Magnetic potential
21 = External electric current 22 = Reaction electric current 23 = Pore pressure
24 = External mass flux 25 = Reaction mass flux 26 = Bearing pressure
27 = Bearing force 28 = Velocity 29 = Rotational velocity
30 = Acceleration 31 = Rotational acceleration 32 = Modal mass
33 = Rotational modal mass 34 = Contact normal stress 35 = Contact normal force
36 = Contact friction stress 37 = Contact friction force 38 = Contact status
39 = Contact touched body 40 = Herrmann variable 41 = Pyrolyzed mass density
42 = Mass rate of gas 43 = Solid density rate 44 = Liquid density rate
45 = Coke density rate 46 = Tying force 47 = Coulomb force
48 = Tying moment

1.5. User subroutine: Examples

Several examples for the usage of user subroutines in Simufact Forming are described in this section. The simulation
models of these examples can be found in the Simufact demos, that can be open through the menu Help > Demos &
examples in Simufact Forming. The used subroutines are placed in the folder AdditionalFiles of each project.

1.5.1. State and post variables

Simufact Forming allows the user to define additional memory for every node and element to store additional infor-
mation during the simulation. This feature becomes useful when, in addition to usual output variables, additional in-
formation of nodal and/or elemental positions is required. One such example could be a total strain in the system.
Simufact Forming by default outputs the elastic and plastic strains separately. Through elemental state variables, the
two can be combined and displayed. These variables are then calculated through a user-defined formula through a
user subroutine.

State variables can be displayed in the result window through activating post variables in the forming control dialog
window. In this example, two nodal and one elemental post variables are defined to be displayed.

32
2023.2 User subroutines State and post variables

Figure 1.11. Activating nodal and elemental post variables

Custom names can be assigned to these variables by clicking on the "edit/show" icon ( ) on the right. The following
window will open where the number of variables and its names can be assigned. These names will be used in the
result tree (under User-defined).

Figure 1.12. User-defined nodal post variable name definition

33
2023.2 User subroutines State and post variables

Figure 1.13. User-defined elemental post variable name definition

The total number of user-defined state variables has also to be specified for the solver. This can be done in the forming
control at Forming control > Advanced > User-defined as shown below:

Figure 1.14. Nodal and element state variables selection

The number of selected element state variables is two even though only one is required. By default the first element
state variable is always assigned to the temperature field. Therefore the first elemental state variable that is available
for the user is from position 2. There is no such reservation in the nodal state variable list. Hence users can use nodal
state variables from position 1.

1.5.1.1. Example
A simple 2D upsetting model is made in Simufact Forming. Two nodal and one elemental user variables are considered.
In this example, the first node variable "user-nodevar" stores the increment number and the second node variable
"outline-variable" stores 10x increment number only at the border nodes. The element variable "Total-strain" stores

34
2023.2 User subroutines State and post variables

the total strain (elastic strain + plastic strain) for every element in the model. These custom variables are listed in the
result window under User-defined.

Figure 1.15. User-defined element and node variables in the result tree
The used user subroutine is shown below. The subroutine is self-explanatory with the written comments. Four sub-
routines are used within it: ueloop() and plotv() to calculate and plot the element variables whereas ufninc() and
upstno() are used for the node variables.

subroutine ueloop(m1,m2,iflag)
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
c ** Start of generated type statements **
#include "hards.cmn" !nstats
#include "creeps.cmn" !timinc
#include "arrays.cmn" !iintel
#include "elemdata.cmn" !ieltype
#include "space.cmn" !ints
#include "concom.cmn" !ipass
c
integer iflag, m1, m2
c ** End of generated type statements **
integer layer, icode,ioflag,ibodyid,nn1,inn,i1,ityp,
* istat
real*8 eplas,deplas
c
c user routine that gets called in key element loops
c may be used in a variety of ways
c
c m1 user element number
c m2 internal element number
c iflag 1 - opress loop
c 2 - oasemb loop
c 3 - oasmas loop
c 4 - ogetst loop
c 5 - scimp loop (end of increment)
c
c IPASS = 1 STRESS
c IPASS = 2 HEAT TRANSFER
c IPASS = 3 FLUIDS
c IPASS = 4 JOULE HEATING
c IPASS = 5 DIFFUSION
c IPASS = 6 ELECTROSTATICS
c IPASS = 7 MAGNETOSTATIC
c IPASS = 8 ELECTROMAGNETICS

if(ipass.ne.1) goto 999 !only for stress pass

if(iflag.eq.5.and.nstats.gt.1)then
c utility routine to extract contact body ID

35
2023.2 User subroutines State and post variables

c
c ielem element number
c ibodyid body ID
c ioflag =0, user element number
c =1, element number is internal
ioflag=0
call get_bodyid(m1,ibodyid,ioflag)
if(ibodyid.ne.1)goto 999
c
c get number of integrationpoints for element (nn1)
ityp=ieltype(m2)
i1=ityp-1
nn1=ints(i1+iintel)
c
do inn=1,nn1
layer=1
icode=28 !strain rate
call get_elmvar(icode,m1,m2,inn,layer,deplas)
c calculate incremental plastic strain
deplas=deplas*timinc
c
c read data ioflag = 0
c write data ioflag = 1
c
c get 2nd state variable (ioflag=0)
istat=2
icode=1000000+istat
ioflag=0
call elmvar1(icode,m1,inn,layer,eplas,ioflag)
c cummulate value
eplas=eplas+deplas
c store 2nd state variable (ioflag=1)
istat=2
icode=1000000+istat
ioflag=1
call elmvar1(icode,m1,inn,layer,eplas,ioflag)
enddo
endif

999 return
end

subroutine plotv(v,s,sp,etot,eplas,ecreep,t,m,nn,kcus,ndi,
* nshear,jpltcd)
c* * * * * *
c
c define a variable for contour plotting (user subroutine).
c
c v variable to be put onto the post file
c s (idss) stress array
c sp stresses in preferred direction
c etot total strain (generalized)
c eplas total plastic strain
c ecreep total creep strain
c t array of state variable (temperature first)
c m(1) user element number
c m(2) internal element number

36
2023.2 User subroutines State and post variables

c m(3) material id
c m(4) internal material id
c nn integration point number
c kcus(1) layer number
c kcus(2) internal layer number
c ndi number of direct stress components
c nshear number of shear stress components
c jpltcd the absolute value of the user's entered post code
c
c
c the components of s, sp, etot, eplas and ecreep are given in the
order
c as defined in volume B for the current element.
c
c* * * * * *
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif

c ** Start of generated type statements **

real*8 ecreep, eplas, etot
integer jpltcd, kcus, m, ndi, nn, nshear
real*8 s, sp, t, v
c ** End of generated type statements **
dimension s(*),etot(*),eplas(*),ecreep(*),sp(*),m(4),kcus(2),
* t(*)
integer ioflag,icode

if(jpltcd.eq.1)then
c get elastic strain
icode=127 !equivalent elastic strain post code
ioflag=0 !read
call elmvar1(icode,m(1),nn,kcus(1),v,ioflag)
v=v+t(2) !add 2nd saved state variable
endif
return
end

c #########################################################

subroutine ufninc(inc,timinc,cptim,ipost,jxtrac)
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
#include "spaceco.cmn"
#include "spacevec.cmn" !inocon

integer inc,ipost,jxtrac
real*8 timinc,cptim
c
c inc increment number
c timinc time step
c cptim time at begining of increment

37
2023.2 User subroutines State and post variables

c ipost current setting to control writing of post flag

c 0 - don't write
c 1 - write
c jxtrac control putting mesh on post file
c if set to 1, then lm(1)=3

integer ind,istat,ioflag,iind,ibd
real*8 var

if(inc.gt.0)then
ibd=1
ioflag=1
c read data ioflag = 0
c write data ioflag = 1
do ind=1,numnp_cn
istat=1
var=real(inc)
c saving "inc" in first node variable
call ndstatvar(istat,ind,var,ioflag)
c---------------------------------------------
istat=2
c determine boundary node, i.e. return value of inocon>0
iind=inocon(ind)
if(iind.gt.0)then
var=real(inc)*10.d0
else
var=0.d0
endif
c saving "10 x inc" in only case of border node variable
call ndstatvar(istat,ind,var,ioflag)
enddo
endif

return
end

subroutine upstno(nqcode,nodeid,valno,nqncomp,nqtype,
* nqaver,nqcomptype,nqdatatype,
* nqcompname)
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
c ** Start of generated type statements **
integer nodeid, nqaver, nqcode, nqcomptype, nqdatatype, nqncomp
integer nqtype
real*8 valno
c ** End of generated type statements **
c
dimension valno(*)
character*24 nqcompname(*)
c
c input: nqcode user nodal post code , e.g. -1
c nodeid node id

38
2023.2 User subroutines State and post variables

c nqcompname not used (future expansion)

c
c output: valno() nodal values
c real/imag valno( 1: nqncomp) real
c valno(nqncomp+1:2*nqncomp) imag
c magn/phas valno( 1: nqncomp) magn
c valno(nqncomp+1:2*nqncomp) phas
c nqncomp number of values in valno
c nqtype 0 = scalar
c 1 = vector
c nqaver only for DDM 0 = sum over domains
c 1 = average over domains
c nqcomptype 0 = global coordinate system (x,y,z)
c 1 = shell (top,bottom,middle)
c 2 = order (1,2,3)
c nqdatatype 0 = default
c 1 = modal
c 2 = buckle
c 3 = harmonic real
c 4 = harmonic real/imaginary
c 5 = harmonic magnitude/phase
c
c
c to obtain nodal values to be used in this subroutine from
c the marc data base the general subroutine NODVAR is available:
c
c call nodvar(icod,nodeid,valno,nqncomp,nqdatatype)
c
c output: valno
c nqncomp
c nqdatatype
c
c input: nodeid
c icod
c 0='Coordinates '
c 1='Displacement '
c 2='Rotation '
c 3='External Force '
c 4='External Moment '
c 5='Reaction Force '
c 6='Reaction Moment '
c 7='Fluid Velocity '
c 8='Fluid Pressure '
c 9='External Fluid Force '
c 10='Reaction Fluid Force '
c 11='Sound Pressure '
c 12='External Sound Source '
c 13='Reaction Sound Source '
c 14='Temperature '
c 15='External Heat Flux '
c 16='Reaction Heat Flux '
c 17='Electric Potential '
c 18='External Electric Charge '
c 19='Reaction Electric Charge '
c 20='Magnetic Potential '
c 21='External Electric Current'
c 22='Reaction Electric Current'
c 23='Pore Pressure '

39
2023.2 User subroutines State and post variables

c 24='External Mass Flux '

c 25='Reaction Mass Flux '
c 26='Bearing Pressure '
c 27='Bearing Force '
c 28='Velocity '
c 29='Rotational Velocity '
c 30='Acceleration '
c 31='Rotational Acceleration '
c 32='Modal Mass '
c 33='Rotational Modal Mass '
c 34='Contact Normal Stress '
c 35='Contact Normal Force '
c 36='Contact Friction Stress '
c 37='Contact Friction Force '
c 38='Contact Status '
c 39='Contact Touched Body '
c 40='Herrmann Variable '
c 41='Density of Solid '
c 42='Mass Flow Rate of Gas '
c 43='Rt of Chnge of Pyrolysis '
c 44='Rt of Chnge of Liquid Den'
c 45='Rt of Chnge of Coke Densi'
c
integer istat,ioflag,nodint,inid
real*8 var

c.....converting external nodeid to internal

c.....For version > 12.0
inid=nodint(nodeid)
if(nqcode.lt.0)then
istat=abs(nqcode)
ioflag=0
c only scalar values are supported so far for ARC
nqtype=0
nqdatatype=0
c read data ioflag = 0
c write data ioflag =1
call ndstatvar(istat,inid,var,ioflag)
valno(1)=var
endif

return
end

This subroutine has to be specified in the forming control as shown in the following figure:

40
2023.2 User subroutines Motion

Figure 1.16. Specify the user subroutine in the forming control

1.5.2. Motion
The user subroutine umotion can be used to define the motion of a body during the simulation. This subroutine should
only be used for bodies with velocity controlled rigid surfaces. This subroutine is called during the calculation at the
beginning of each time increment and returns the surface velocities for that increment. Also refer to Marc Volume
D: User Subroutines and Special Routines - User-defined Loading, Boundary Conditions, and State Variables User
Subroutines List for more details.

1.5.2.1. Example
A simple 2D and 3D upsetting model is set up in Simufact Forming. With the help of the subroutine umotion(),
the deformation process will be stopped when the force on the top die exceeds the user-defined force of 6kN. The
movement of the top die will then be reversed with the user-defined velocity of 20mm/sec. The control information,
in this case die force and die velocity, are entered through the user-defined variables defined in the forming control
under Advanced > User-defined. The same subroutine is used for both simulation types (2D and 3D). Additionally
the box for the subroutine umotion has to be checked in the tab User subroutine.

41
2023.2 User subroutines Motion

Figure 1.17. User-defined real variables and subroutine in the forming control

Figure 1.18. Select the subroutine umotion in the forming control

The user subroutine for reversing the forming process in case that the force exceeds the user-defined value is given
below:

subroutine motion(x,f,v,time,dtime,nsurf,inc)
c
c
c * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
c
c user routine to provide surface motion data
c
c------------------------------------------------------------------------
c Caution: The increment number (inc), which is passed from the
c calling routine, is useful in making sure the data given
c is appropriately attributed to the relevant increments.
c For example, a center of rotation given for increment 0
c if(inc.eq.0) then
c x(1)=...
c etc.
c endif
c will be updated internally as motion and deformation take
c place. However, if the increment is not specified, then
c this center of rotation will be reset each time this routine
c is read, i.e. every increment. Obviously the results will
c be different for the two cases.
c
c The same caution is also applicable to the surface number
c (nsurf). Thus the data should be coupled only with the
c surface number(s) for which it is meant. For example:

42
2023.2 User subroutines Motion

c if(nsurf.eq.2 .or. nsurf.eq.4) then

c x(1)=...
c ...
c v(1)=...
c etc.
c endif
c
c Important: The correct use of increment numbers and surface numbers
c may give rise to more complex (or nested) if statements.
c For example:
c if(inc.eq.2 .and. nsurf.eq.3) then
c ...
c elseif(....) then
c if(....) then
c ...
c endif
c ...
c endif
c------------------------------------------------------------------------
c
c 2-d:
c input : nsurf - number of the surface for which data is
c requested
c time - the time at which data is requested
c dtime - the current time increment
c x(3) - current die defining coordinates:
c x(1) = 1st coordinate of center of
c rotation
c x(2) = 2nd coordinate of center of
c rotation
c x(3) = angle rotated around z-axis
c
c f(3) - the current surface load:
c f(1) = 1st component of load
c f(2) = 2nd component of load
c f(3) = moment
c inc - the increment number
c
c output : v(3) - current surface velocities
c v(1) = 1st component of the velocity at
c the center of rotation.
c v(2) = 2nd component of the velocity at
c the center of rotation.
c v(3) = angular velocity
c
c 3-d:
c input : nsurf - the number of the surface for which data
c is requested
c time - the time at which data is requested
c dtime - the current time increment
c x(6) - current die defining coordinates:
c x(1) = 1st coordinate of center of
c rotation
c x(2) = 2nd coordinate of center of
c rotation
c x(3) = 3rd coordinate of center of
c rotation
c Axis for specifying angular velocity:

43
2023.2 User subroutines Motion

c x(4) = 1st component of direction cosine

c x(5) = 2nd component of direction cosine
c x(6) = 3rd component of direction cosine
c
c f(6) - the current surface load:
c f(1) = 1st component of load
c f(2) = 2nd component of load
c f(3) = 3nd component of load
c f(4) = 1st component of moment
c f(5) = 2nd component of moment
c f(6) = 3rd component of moment
c inc - the increment number
c
c output : v(4) - current surface velocities
c v(1) = 1st component of the velocity at
c the center of rotation
c v(2) = 2nd component of the velocity at
c the center of rotation
c v(3) = 3nd component of the velocity at
c the center of rotation
c v(4) = angular velocity around axis defined
c above with x(4), x(5), and x(6).
c
c
c * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
c
c
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
#include "spaceco.cmn" !contact body information, cbname, numdie_cn
#include "machin.cmn" !iprtl
#include "prepro.cmn" !kou
#include "sf_userdata.cmn"
c ** Start of generated type statements **
real*8 dtime, f
integer inc, nsurf
real*8 time, v, x
c ** End of generated type statements **
dimension x(*),v(*),f(*)
character*24 cbody,cbody0
integer i,imvbody,i2switch,ilen1
real*8 mxforce,vel
c
c
c first find out with contact body ID should be moved
if(inc.eq.0.and.nsurf.eq.1)then
i2switch=0
cbody='Punch'
ilen1=index(cbody,' ') - 1
c nruser - number of real user variables
if(nruser.lt.2)then
write(iprtl,1001) nruser,2
call quit(13)
endif
c get information from user data block check units!!

44
2023.2 User subroutines Damage

c in forming SI-m unit are used

mxforce=ruserdata(2)*1000.d0
vel=ruserdata(1)/1000.d0
c get body id
do i=1,numdie_cn
cbody0=cbname(i)
if(cbody0(:ilen1).eq.cbody(:ilen1))then
imvbody=i
exit
endif
enddo
endif

c skip alll bodies which should not move

if(nsurf.ne.imvbody)goto 999
c
c check if max force is reached
if(i2switch.eq.0)then
if(abs(f(1)).gt.mxforce)then
write(kou,1002)mxforce, cbody(:ilen1),abs(f(1))
i2switch=1
endif
endif
c define velocity
v(1)=vel
if(i2switch.eq.1)v(1)=-v(1)
write(kou,1003)cbody(:ilen1),v(1)
999 return
1001 format(13x,'*** error ',
* 'Number of user defined real variables nruser:',i3,'
* is less then needed: ',i3,/13x,'check input deck!')
1002 format(13x,'max force of ',f12.6,' is reached in body:'a,
* ' force is:',f12.6)
1003 format(13x,'Calculated velocity of Body ',a,' is:',f12.6)

end

1.5.3. Damage
Damage information during the simulation can be used to predict the failure of the workpiece. This is done by com-
paring the calculated solution to some failure criteria. Simufact Forming by default supports seven damage models,
i.e. Cockroft-Latham, Modified Mohr-Coulomb, Lemaitre, Oyane, Johnson-Cook, Gurson and Bonora.

Figure 1.19. Supported damage models in Simufact Forming

However only one of them can be selected as output in the simulation. With the help of a user subroutine, multiple
user-defined damage models can be implemented simultaneously in the simulation.

45
2023.2 User subroutines Damage

1.5.3.1. Example
In this example, damage values through three different models are calculated and stored in three elemental state vari-
ables. In addition to linking the subroutine, three more elemental state variables and three user-defined elemental post
variables have to be selected and defined in the forming control, as shown in the figures below:

Figure 1.20. Subroutine, state and post variable definition

Figure 1.21. Definition of element post variables

The example consist of two stages over stage control that share the same user subroutine, as written below:

subroutine ueloop(m1,m2,iflag)
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)

46
2023.2 User subroutines Damage

#endif

c* * * * * *
c..................
c Output values for the integral damage criteria
c damage growing at positive stress only
c..................
c 2 [Atk81] mean normal stress
c A. G. Atkins, 1981, Possible explanation for unexpected
c departures in hydrostatic tension–fracture strain
c relations. Metal Science, pp. 81–83
c..................
c 3 [AyHi87] triaxiality
c M. Ayada, T. Higashino and K. Mori, 1987,
c Central Bursting in Extrusion of Inhomogeneous Materials.
c Advanced Technology of Plasticity, Vol. 1, pp. 553-558
c..................
c 4 [OySa80] modified triax
c M. Oyane, T. Sato, K. Okimoto and S. Shima, 1980,
c Criteria for ductile fracture and their applications,
c J. Mech. Work. Technol, Vol. 4, pp. 65-81
c..................
c* * * * * *

c ** Start of generated type statements **

#include "hards.cmn" !nstats
#include "space.cmn" !ints
#include "creeps.cmn" !timinc
#include "concom.cmn" !ipass
#include "prepro.cmn" !kou
#include "machin.cmn" !iprtl
#include "arrays.cmn" !iintel
#include "elemdata.cmn" !ieltype
#include "marc_prestt.cmn" !iprestt_data

integer iflag, m1, m2

c ** End of generated type statements **

integer layer, icode,ioflag,ibodyid,nn1,inn,i1,ityp,

* istat,ndi,nshear
real*8 eplas,deplas,val,vmises,hydro,triax,erate
real*8 depsilon,ddamage,biggg,c
c
c user routine that gets called in key element loops
c may be used in a variety of ways
c
c m1 user element number
c m2 internal element number
c iflag 1 - opress loop
c 2 - oasemb loop
c 3 - oasmas loop
c 4 - ogetst loop
c 5 - scimp loop (end of increment)

c IPASS = 1 STRESS
c IPASS = 2 HEAT TRANSFER
c IPASS = 3 FLUIDS
c IPASS = 4 JOULE HEATING

47
2023.2 User subroutines Damage

c IPASS = 5 DIFFUSION
c IPASS = 6 ELECTROSTATICS
c IPASS = 7 MAGNETOSTATIC
c IPASS = 8 ELECTROMAGNETICS

if(ipass.ne.1) goto 999 !only for stress pass

if(inc.eq.0)then
c check first if total number is set correct for the job
if(nstats.lt.4)then
write(kou,1002)nstats,4
write(iprtl,1002)nstats,4
call quit(13)
endif

c initialization of the state variables

if(iprestt_data(1).ne.0) goto 100
c get number of integrationpoints for element (nn1)
ityp=ieltype(m2)
i1=ityp-1
nn1=ints(i1+iintel)
do inn=1,nn1
ioflag=1 !write data
layer=1
do istat=2,4
icode=1000000+istat
c initialize to zero only if iprestt_data(1).eq.0 (no multistage job)
c not for a multistage job
val=0.0d0 !initialize to zero
call elmvar1(icode,m1,inn,layer,val,ioflag)
enddo
enddo
goto 999
endif

100 if(inc.eq.0) goto 999

if(iflag.eq.5.and.nstats.gt.1)then

c utility routine to extract contact body ID

c ielem element number
c ibodyid body ID
c ioflag =0, user element number
c =1, element number is internal
ioflag=0
call get_bodyid(m1,ibodyid,ioflag)
if(ibodyid.ne.1)goto 999
c
c get number of integrationpoints for element (nn1)
ityp=ieltype(m2)
i1=ityp-1
nn1=ints(i1+iintel)

do inn=1,nn1
layer=1
icode=17 !vonmises stress
call get_elmvar(icode,m1,m2,inn,layer,vmises)

48
2023.2 User subroutines Damage

icode=18 !hydrostatic stress

call get_elmvar(icode,m1,m2,inn,layer,hydro)

triax=hydro/vmises

if(triax.ge.0) then

c equivalent plastic strain increment delta_epsilon

c delta_epsilon = (d_epsilon/d_t) * delta_t
icode = 28 !plastic strain
call get_elmvar(icode,m1,m2,inn,layer,erate)
depsilon=erate*timinc

c.......................................................................
c2 Damage criteria by [Atk81] 1
istat=2
icode=1000000+istat
ioflag=0 !read previous value
call elmvar1(icode,m1,inn,layer,val,ioflag)

ddamage = hydro * depsilon

val=val+ ddamage

ioflag=1 !write
call elmvar1(icode,m1,inn,layer,val,ioflag)

c.......................................................................
c3 Damage criteria by [AyaH87]
istat=3
icode=1000000+istat
ioflag=0 !read previous value
call elmvar1(icode,m1,inn,layer,val,ioflag)

ddamage = (hydro/vmises) * depsilon

val=val+ ddamage

ioflag=1 !write
call elmvar1(icode,m1,inn,layer,val,ioflag)

c.......................................................................
c4 Damage criteria by [OySa80]
istat=4
icode=1000000+istat
ioflag=0 !read previous value
call elmvar1(icode,m1,inn,layer,val,ioflag)

c=0.424
ddamage=(1+(hydro/(c*vmises)))*depsilon
val=val+ ddamage

ioflag=1 !write
call elmvar1(icode,m1,inn,layer,val,ioflag)
c.......................................................................
endif
enddo
endif

49
2023.2 User subroutines Damage

999 return
1002 format(13x,'*** error: Total number of state variables ',i5,
* ' is less than needed in plotv:',i5,
* /13x,'increase number to "element state variables" in
* job options under:',
* /13x,'Forming control-->Advanced-->Solver-->Advanced ',/)

end

subroutine plotv(v,s,sp,etot,eplas,ecreep,t,m,nn,kcus,ndi,
* nshear,jpltcd)
c* * * * * *
c
c define a variable for contour plotting (user subroutine).
c
c v variable to be put onto the post file
c s (idss) stress array
c sp stresses in preferred direction
c etot total strain (generalized)
c eplas total plastic strain
c ecreep total creep strain
c t array of state variable (temperature first)
c m(1) user element number
c m(2) internal element number
c m(3) material id
c m(4) internal material id
c nn integration point number
c kcus(1) layer number
c kcus(2) internal layer number
c ndi number of direct stress components
c nshear number of shear stress components
c jpltcd the absolute value of the user's entered post code
c
c
c the components of s, sp, etot, eplas and ecreep are given in the
order
c as defined in volume B for the current element.
c
c* * * * * *
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
c ** Start of generated type statements **
real*8 ecreep, eplas, etot
integer jpltcd, kcus, m, ndi, nn, nshear
real*8 s, sp, t, v
c ** End of generated type statements **
dimension s(*),etot(*),eplas(*),ecreep(*),sp(*),m(4),kcus(2),
* t(*)

if(jpltcd.eq.2) v=t(2) !1st user state variable

if(jpltcd.eq.3) v=t(3) !2nd user state variable
if(jpltcd.eq.4) v=t(4) !3rd user state variable

return

50
2023.2 User subroutines Friction

end

In the results some negative damage values may occur. During remeshing the results of the integration points of the
old elements have to be transferred to the integration points of the new elements. This is done by interpolation and
if needed geometrically by extrapolation. This leads to some unavoidable inaccuracy and in this case to the negative
damage values.

The user may add a check to the user subroutine to force the calculated damage values to be equal or bigger than zero.

1.5.4. Friction
Through the subroutine ufric, the friction behavior of the contact surfaces between the die and the workpiece can be
changed dynamically. The ufric subroutine is called when either only shear or coulomb friction is added in the process
tree. The subroutine is called for outer elements and all nodes that lie on the workpiece surface. The subroutine only
works for node to segment contact.

1.5.4.1. Example
In the 2D upsetting example, the starting roughness value of the material is input through a real user variable in the
forming control. Based on the roughness value of the surface node, the friction coefficient is chosen in the ufric
subroutine. At the end of each increment, the roughness values for the surface nodes are changed as a function of
von Mises and contact normal stress. The changed roughness value is then used to set the friction coefficient in the
simulation for the next increment. For using ufric subroutine, an extra subroutine flag needs to be activated in the
GUI. The definition of the user subroutine and the nodal state and post variables are shown in the figures below:

Figure 1.22. Enabling ufric subroutine in Forming control

51
2023.2 User subroutines Friction

Figure 1.23. Subroutine, state and post variable definition

Figure 1.24. Definition of nodal post variables

The subroutine is shown below:

subroutine ufric(mibody,x,fn,vrel,temp,yiel,fric,time,inc,nsurf)
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)

52
2023.2 User subroutines Friction

#endif

#include "prepro.cmn" ! kou

#include "spaceco.cmn" ! jfnold
#include "sf_userdata.cmn" ! ruserdata

c ** Start of generated type statements **

real*8 x, fn, fric
integer inc, mibody, nsurf
real*8 temp, time, vrel, yiel
c ** End of generated type statements **
dimension mibody(*),x(*),vrel(*),temp(*)
c
c user routine to define friction behavior
c
c Input
c for distributed friction based upon nodal stress
c
c mibody(1) - user element number
c mibody(2) - edge or face number
c mibody(3) - edge/face integration point number
c mibody(4) - internal element number
c fn - normal stress
c
c for nodal friction based upon nodal forces
c
c mibody(1) - user node number
c mibody(2) - not used
c mibody(3) - not used
c mibody(4) - internal node number
c fn - normal force
c
c for either model
c
c x - updated coordinates of contact point where
c friction is being calculated
c vrel - relative sliding velocity in tangential direction(s)
c temp(1) - temperature
c temp(2) - voltage in Joule heating simulation
c yiel - flow stress
c time - time at begining of increment
c inc - increment number
c nsurf - surface being contacted by the side for
c which friction calculation is being made
c
c Output
c fric - friction coefficint or factor
c

c skip symmetry body

c jfnold(4,i) = 1 : rigid
c 2 : deformable
c 3 : symmetry
c 4 : heat-rigid
c 5 : workpiece (autoforge only)

real*8 xrough,xrough0
integer istat,ioflag

53
2023.2 User subroutines Friction

if(jfnold(4,nsurf).eq.3)goto 999

c read data ioflag = 0, write data ioflag = 1

ioflag=0
istat=1
call ndstatvar(istat,mibody(4),xrough,ioflag)
if(xrough.gt.0)then
if(xrough.lt.30.0) then
fric=0.15
elseif (xrough.lt.35.0) then
fric=0.25
elseif (xrough.lt.50.0) then
fric=0.40
else
fric=0.50
end if
endif
999 return
end

c #############################################################
subroutine uedinc(inc,incsub)
c called from oscinc

#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif

#include "prepro.cmn" ! kou

#include "machin.cmn" ! iprtl
#include "spaceco.cmn" ! jfnold
#include "spacevec.cmn" ! inocon
#include "sf_userdata.cmn" ! ruserdata

c
c... subroutine which is called at the end of each
c... increment
c...
c... inc : increment number
c... incsub: sub-increment number
c
c
c-------------------------------------------------------
c subroutine at the end of each inc
c can be used to initialize nodal state variables at inc 0
c-------------------------------------------------------
c ** Start of generated type statements **
integer inc, incsub
c ** End of generated type statements **
integer ind,istat,ioflag,iind
real*8 xrough,var

xrough=0.0d0 ! initialize
if(inc.eq.0)then

54
2023.2 User subroutines Friction

c
c check first if nodal state variable is defined
c nstatev = number of nodal state variables
c nruser - number of real user variables
c ruserdata(nruser) - real user data array
c
if(nstatev.gt.0)then
if(nruser.gt.0)then
xrough=ruserdata(1)
else
write(kou,1002)nruser,1,xrough
call quit(13)
endif
else
write(iprtl,1001)nstatev,1
write(kou,1001)nstatev,1
call quit(13)
endif

do ind=1,numnp_cn
c figure out if node is a boundary node
c do nn=1,numnp
c inn=inocon(nn)
c if(inn.lt.0) then
c check if Internal Node belongs to body ibd only 2D
c if(abs(inn).eq.ibd) then
c endif
c elseif (inn.gt.0) then
c inn1=mod(inn,nbmx)
cc check if Boundary Node belongs to body ibd
c if (inn1.eq.ibd) then
c endif
c endif

iind=inocon(ind)
if(iind.gt.0)then
var=xrough
else
var=0.d0
endif

c read data ioflag = 0, write data ioflag = 1

ioflag=1
istat=1
call ndstatvar(istat,ind,var,ioflag)
enddo
endif !end inc loop
999 return

1001 format(13x,'*** error ',

* 'Number of user defined nodal variables nstatev:',i3,'
* is less then needed: ',i3,/13x,'check input deck!')
1002 format(13x,'*** error ',
* 'Number of user defined real variables nruser:',i3,'
* is less then needed: ',i3,/13x,'fixed value is used: ',
* f12.6,' check input deck!')
end

55
2023.2 User subroutines Friction

subroutine ufninc(inc,timinc,cptim,ipost,jxtrac)
c
c ###########################################################
c usersubroutine at the end of an increment before result output.
c ###########################################################
c
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif

#include "dimen.cmn" !ncrd

#include "prepro.cmn" ! kou
#include "spaceco.cmn" ! jfnold
#include "spacevec.cmn" ! area_d
#include "sf_userdata.cmn" ! ruserdata

c ** Start of generated type statements **

integer inc,ipost,jxtrac
real*8 timinc,cptim
c ** End of generated type statements **
c
c inc increment number
c timinc time step
c cptim time at begining of increment
c ipost current setting to control writing of post flag
c 0 - don't write
c 1 - write
c jxtrac control putting mesh on post file
c if set to 1, then lm(1)=3
integer ind,iflag,istat,lnode,lnpost,lepost,num_node,
* num_npost,num_epost,nsize_node,ioflag
real*8 xrough,xrough0,vmises,cstress,var_node,var_elem,x1
real*8 contact_area,cnforce
dimension lnode(1),lnpost(1),lepost(1),var_node(3),var_elem(1)
c ----------------------------------------------------------------------
c subroutine get_nodvar(lnode,lnpost,lepost,var_node,var_elem,
c & num_node,num_npost,num_epost,nsize_node,
c & iflag)
c
c routine to extract nodal values based on the user node numbers,
c nodal post codes and/or element post codes.
c for element post codes, the integration values will be extrapolated to
c the nodes and averaged.
c lnode() - user node numbers
c lnpost() - nodal post codes
c lepost() - element post codes

c var_node() - return nodal results from nodal post codes

c var_elem() - return nodal results from element post codes
c num_node - number of nodes
c num_npost - number of nodal post codes
c num_epost - number of element post codes
c nsize_node - dimension of var_node()
c For example, for a post-code “displacement”
c you will get 3 nodal displacements in X,Y,Z.

56
2023.2 User subroutines Friction

c In this case, num_npost=1 and nsize_node=3.

c iflag - =0 node numbers are external
c =1 node numbers are internal
c----------------------------------------------------------------------------------

c skip calculation for increment 0

if(inc.eq.0) goto 999
xrough0=ruserdata(1)
xrough=0.0d0
do ind=1,numnp_cn
c -- read data ioflag = 0, write data ioflag = 1
ioflag=0
istat=1
call ndstatvar(istat,ind,xrough,ioflag) !get stored roughness
if(xrough.gt.0)then !then boundary node which roughness
lnode(1)=ind
lnpost(1)=35 !35 Contact Normal Force
lepost(1)=17 !17 Equivalent von Mises stress
num_node=1
num_npost=1
num_epost=1
nsize_node=3
iflag=1 !node numbers are internal
call get_nodvar(lnode,lnpost,lepost,var_node,
& var_elem,num_node,num_npost,num_epost,nsize_node,
& iflag)
vmises=var_elem(1)

c -- Calculating contact normal stress

contact_area=area_d(ind) !get nodal area
if(contact_area.ne.1.0d+20) then
cnforce=var_node(1)**2+var_node(2)**2
if(ncrd.eq.3) cnforce = cnforce + var_node(3)**2
cnforce=sqrt(cnforce)
cstress=cnforce/contact_area
if(vmises.gt.0.0d0 .and. cstress.ne.0.0d0) then
x1=cstress/vmises
x1=abs(xrough0*(1 - x1))

if(x1.lt.xrough) then !if roughness is smaller save it

ioflag=1
istat=1
call ndstatvar(istat,ind,x1,ioflag)
endif
endif
endif
endif
enddo
999 return
end

c ###########################################################
c subroutine for result output
c ###########################################################
subroutine upstno(nqcode,nodeid,valno,nqncomp,nqtype,
* nqaver,nqcomptype,nqdatatype,
* nqcompname)

57
2023.2 User subroutines Friction

#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif

c ** Start of generated type statements **

integer nodeid, nqaver, nqcode, nqcomptype, nqdatatype, nqncomp
integer nqtype
real*8 valno
c ** End of generated type statements **
c
dimension valno(*)
character*24 nqcompname(*)
c
c input: nqcode user nodal post code , e.g. -1
c nodeid node id
c nqcompname not used (future expansion)
c
c output: valno() nodal values
c real/imag valno( 1: nqncomp) real
c valno(nqncomp+1:2*nqncomp) imag
c magn/phas valno( 1: nqncomp) magn
c valno(nqncomp+1:2*nqncomp) phas
c nqncomp number of values in valno
c nqtype 0 = scalar
c 1 = vector
c nqaver only for DDM 0 = sum over domains
c 1 = average over domains
c nqcomptype 0 = global coordinate system (x,y,z)
c 1 = shell (top,bottom,middle)
c 2 = order (1,2,3)
c nqdatatype 0 = default
c 1 = modal
c 2 = buckle
c 3 = harmonic real
c 4 = harmonic real/imaginary
c 5 = harmonic magnitude/phase
c
c
integer istat,ioflag,nodint,inid
c
c.....converting external nodeid to internal
c.....For version 12.0 and higher
inid=nodint(nodeid)
valno(1)=0.d0
nqtype=0
nqdatatype=0
nqncomp=1

if(nqcode.eq.-1)then !first user post variable

c read data ioflag = 0, write data ioflag = 1
ioflag=0
istat=1
call ndstatvar(istat,inid,valno(1),ioflag)
endif
return
end

58
2023.2 User subroutines Table evaluation

1.5.5. Table evaluation

In this example, the utility routine tabva2() is used to calculate the values from the input data table used in the sim-
ulation through user subroutines. Every input information provided in a table should have a unique table ID for the
solver to distinguish between different tables. The tabva2() utility routine helps the user to calculate information from
the table within the subroutine. Users can use the evaluated value further in their calculation.

1.5.5.1. Example
A simple upsetting job is set up. To know the table id given by the Simufact Forming GUI, one should temporarily
start and stop the simulation and browse for the job *.dat- and material *.umt-file created in the _Run_ directory of
the process. Open the material input file (*.umt) in a text editor to see that various material parameters are provided in
a table as a function of the temperature. The input information about the Young's modulus from the example material
file (*.umt) is shown below:

$---------------------------------------------------------------
$ youngs_modulus
$---------------------------------------------------------------
table, youngs_modulus
$ table id, nb of independent variables, unit number to read data, disable
output, method to read function
11000,1,0,1,0,
$ temperature, nb of points, deny extrapolation
12,24,1,
$ temperature data
3.2315000000000+02,3.7315000000000+02,4.2315000000000+02,4.7315000000000+02,
5.2315000000000+02,5.7315000000000+02,6.2315000000000+02,6.7315000000000+02,
7.2315000000000+02,7.7315000000000+02,8.2315000000000+02,8.7315000000000+02,
9.2315000000000+02,9.7315000000000+02,1.0231500000000+03,1.0731500000000+03,
1.1231500000000+03,1.1731500000000+03,1.2231500000000+03,1.2731500000000+03,
1.3231500000000+03,1.3731500000000+03,1.4231500000000+03,1.4731500000000+03,
$ youngs_modulus data
2.1000000000000+11,2.0700000000000+11,2.0300000000000+11,2.0000000000000+11,
1.9600000000000+11,1.9300000000000+11,1.8800000000000+11,1.8400000000000+11,
1.7900000000000+11,1.7500000000000+11,1.6900000000000+11,1.6400000000000+11,
1.5800000000000+11,1.5300000000000+11,1.4700000000000+11,1.4000000000000+11,
1.3400000000000+11,1.2700000000000+11,1.2000000000000+11,1.1300000000000+11,
1.0600000000000+11,9.8000000000000+10,9.0000000000000+10,8.2000000000000+10,

In this example, we consider the following three tables:

• table id = 11000 : Young's modulus

• table id = 11001 : Thermal expansion coefficient

• table id = 11002 : Thermal conductivity

With the help of the subroutine the current values of these temperature dependent material data will be calculated and
reported at every increment during the simulation. Moreover the Young's modulus value will be converted from the
unit Pascal to Psi. The following settings are made in Forming control > Advanced > User-defined:

• Including the user subroutine and integer user variables

59
2023.2 User subroutines Table evaluation

Figure 1.25. Table id as input through integer user variables

• Using real user variables as scaling factor

Figure 1.26. Scaling factor as input through real user variables

• Check the user-defined element post variables and select three element variables and rename them accordingly

Figure 1.27. Check elemental post variables

60
2023.2 User subroutines Table evaluation

Figure 1.28. Elemental post variables selection

The following subroutine was used to calculate the three table values and display the output. For this example only
the plotv subroutine is used. A user can implement the tabva2 utility in most of the user subroutines. The common
file, ctable.cmn, contains the list of all the independent variables that could be used for the table. In this example,
temperature is used as an independent variable in input hence the variable tempi is updated with the temperature for
the solver to extract the corresponding value.

c--------------------------------------------------------------
c Subroutine to show uses of tabva2 utility routine to extract
c user defined table information provided in the input during
c simulation
c--------------------------------------------------------------
subroutine plotv(v,s,sp,etot,eplas,ecreep,t,m,nn,kcus,ndi,
* nshear,jpltcd)
c* * * * * *
c define a variable for contour plotting (user subroutine).
c
c v variable to be put onto the post file
c s (idss) stress array
c sp stresses in preferred direction
c etot total strain (generalized)
c eplas total plastic strain
c ecreep total creep strain
c t array of state variable (temperature first)
c m(1) user element number
c m(2) internal element number
c m(3) material id
c m(4) internal material id
c nn integration point number
c kcus(1) layer number
c kcus(2) internal layer number
c ndi number of direct stress components
c nshear number of shear stress components
c jpltcd the absolute value of the user's entered post code
c
c
c the components of s, sp, etot, eplas and ecreep are given in the
order
c as defined in volume B for the current element.
c* * * * * *

#ifdef _IMPLICITNONE

61
2023.2 User subroutines Restart with global variables

implicit none
#else
implicit logical (a-z)
#endif

#include "ctable.cmn" !provides variable tempi

#include "sf_userdata.cmn" !provides user data

c ** Start of generated type statements **

real*8 ecreep, eplas, etot
integer jpltcd, kcus, m, ndi, nn, nshear
real*8 s, sp, t, v
c ** End of generated type statements **
dimension s(*),etot(*),eplas(*),ecreep(*),sp(*),m(4),kcus(2),
* t(*)

integer itableid
real*8 orig_tempi, cfactor

c -- check if number of post output is less or equal to the

c -- number of user provided in table id
if(jpltcd.le.niuser) then

c -- Copy user provided table id

itableid=iuserdata(jpltcd)

c -- Copy original tempi to orig_tempi for restoration later

orig_tempi=tempi
tempi=t(1) !First state varibale is the temperature
cfactor=1.0

c -- Selecting multiplier if provided by the user

if(jpltcd.le.nruser) cfactor=ruserdata(jpltcd)

c -- calling the table evaluation routine. The result is saved in

c -- variable v which is put in post file
call tabva2(cfactor,v,itableid,0,0)

c -- Restoring the variable tempi to its original value

tempi=orig_tempi
endif

return
end

1.5.6. Restart with global variables

The recommended way for saving global variables across solver stops is via the
<processname>_kin.res-file introduced in this chapter. The file is handled automatically by the
GUI retaining all expected functionalities of Simufact GUI such as saving the project with results,
Simufact Remote or Restart. In case of a solver stop, the <processname>_kin.res-file is copied from
the current _Run_-directory to the _Results_\Collective-folder.

The solver process might be stopped at any time, e.g. through user intervention, simulation errors or for remeshing
purpose. To properly continue the simulation from a desired restart point all required history dependent information
that is needed in the next steps or increment has to be available and hence must be saved in advance. This also holds for

62
2023.2 User subroutines Restart with global variables

variables which are defined in user subroutines. There are several implemented options to ensure that this is the case.
For example when using state variables (see chapter state variables and post variables) the handling of the accessibility
is ensured by the solver. But when using own common blocks, global variables in modules or even local variables in a
subroutine, this behaviour has to be implemented by the user. Simufact Forming provides an option to store variables
in a file, with the specific ending "kin.res", which is then automatically stored in the results folder when the simulation
aborts. The information therefore remains available for the restart and only hast to be read back accordingly to their
variables.

1.5.6.1. Example
The provided example shows how to share local variables by making them global through a module and additionally
saving the content to *kin.res-file for restart purposes. Making local variables global through common blocks or
modules is a general way for sharing variables among other subroutines or transfer the results throughout multiple
calls of the same subroutine. Another method for sharing is to write the information to a file and read it back, as used
in the motion example earlier. Both methods however get problematic when a restart should be considered:

Global variables in modules are lost when the solver stops. This can already be observed for a normal simulation
when no restart behaviour is required but DDM is used. For a DDM job the simulation is internally stopped as soon as
remeshing occurs and restarted subsequently. The second method - writing information to a general file and reading
it back - will fail for restart as the Simufact Forming GUI will cleanup the run directory at the restart and so the
information gets lost. In both cases the *kin.res-file solution should be applied. The Simufact Forming GUI will
backup the kin.res file before cleaning the directory and will restore them for restart. Also if the job is copied (with
results), the kin.res file is copied, too and is available for possible restart of the job. The usage of the *kin.res-file is
demonstrated in this example.

There are two steps demonstrated: Sharing variables between subroutine calls and supporting of restart. At the end
of each increment the information will be written to the file and will be read at the beginning in case there was a
remeshing. Information that is unique for an increment and is needed for calculations in next increment, irrespective
of nodes/elements, should be saved through this method. If the user needs to store information that is node/element
specific, state variables should be used instead. The handling of state variables in case of a restart is provided by the
solver automatically. The information stored at the state variables will also be mapped to the new mesh by the solver at
the remeshing. Note: In case of result transfer between a multi stage process, the state variables are also automatically
transferred to the next stage, as long as the identical number of state variables is defined in the Forming dialog for
each process.

In this example a simple upsetting model is used. Through the subroutine, increment information is saved in a file
*kin.res. The Simufact Forming GUI will backup and replace the file if there is a job restart. If the job is restarted
from an earlier increment, the correct saved information will be retrieved by the subroutine. Users should modify the
routine depending upon the number of variables to be saved. For demonstration, a real and an integer variable are
assigned and initialized at the beginning. The variables are then linearly increased by factor of one and ten respectively
and plotted in two post variables. Users should stop the running job and restart it from an earlier increment to see if
the information are retained in the restart run. The values of the two post variables can be selected in the result view
under User-defined.

The figures below show the definition of the subroutine and of the post variables in the forming control:

63
2023.2 User subroutines Restart with global variables

Figure 1.29. Subroutine and post variables definition

Figure 1.30. Element post variables definition

The used subroutine is shown below:

c -- Module definition
module uservariables
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
c Define number of real and integer variable to save
integer,parameter :: number_real = 1
integer,parameter :: number_integer = 1

integer,dimension (number_integer) :: integer_array

real,dimension (number_real) :: real_array
end module

64
2023.2 User subroutines Restart with global variables

c -- Read the saved variable at start of increment

subroutine ubginc(inc,incsub)
use uservariables

#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
#include "cdominfo.cmn"
#include "marc_prestt.cmn"

integer inc, incsub

c Initializing variable at simulation start
if(inc.eq.0.and.iprestt_data(23).eq.0) then
real_array(1)=0.0
integer_array(1)=0
endif
call readuservariables !in case of remeshing/Job restart
return
end

c -- Save the variable at end of increment

subroutine ufninc(inc,timinc,cptim,ipost,jxtrac)
use uservariables
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
integer inc,ipost,jxtrac
real*8 timinc,cptim
c inc increment number
c timinc time step
c cptim time at begining of increment
c ipost current setting to control writing of post flag
c 0 - don't write
c 1 - write
c jxtrac control putting mesh on post file
c if set to 1, then lm(1)=3
integer ind,istat,ioflag,iind,ibd
real*8 var

if(inc.gt.0)then
real_array(1)=real_array(1)+1.0
integer_array(1)=integer_array(1)+10
endif
call saveuservariables ! Save the current state to file
return
end

c -- Plot the variables in this example

subroutine plotv(v,s,sp,etot,eplas,ecreep,t,m,nn,kcus,ndi,
* nshear,jpltcd)
use uservariables
c* * * * * *

65
2023.2 User subroutines Restart with global variables

c
c define a variable for contour plotting (user subroutine).
c
c v variable to be put onto the post file
c s (idss) stress array
c sp stresses in preferred direction
c etot total strain (generalized)
c eplas total plastic strain
c ecreep total creep strain
c t array of state variable (temperature first)
c m(1) user element number
c m(2) internal element number
c m(3) material id
c m(4) internal material id
c nn integration point number
c kcus(1) layer number
c kcus(2) internal layer number
c ndi number of direct stress components
c nshear number of shear stress components
c jpltcd the absolute value of the user's entered post code
c
c
c the components of s, sp, etot, eplas and ecreep are given in the
order
c as defined in volume B for the current element.
c
c* * * * * *
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
c ** Start of generated type statements **
real*8 ecreep, eplas, etot
integer jpltcd, kcus, m, ndi, nn, nshear
real*8 s, sp, t, v
c ** End of generated type statements **
dimension s(*),etot(*),eplas(*),ecreep(*),sp(*),m(4)
dimension kcus(2),t(*)

if(jpltcd.eq.1) v=real_array(1)
if(jpltcd.eq.2) v=integer_array(1)
return
end

c -- Function to read saved variables after remeshing

subroutine readuservariables
use uservariables
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
#include "op1.cmn"
#include "sfkin.cmn"
#include "concom.cmn"
#include "cdominfo.cmn"
#include "marc_prestt.cmn"

66
2023.2 User subroutines Restart with global variables

integer i,idone
integer :: iunit = 18
integer :: kin_unit = 57
real*8 dummy
character*(MAX_CARD_DIM_2) filen,filen2

c -- If called for first time then create the file

if(inc.eq.0.and.iprestt_data(23).eq.0.and.iparent.eq.0) then
call get_kin2res_filename(filen, filen2)
call openfile(iunit,filen,1,2)
call sfkin2res (iunit,kin_unit,0,0,0,1,2) !write
call closefile(iunit,filen,0,2)
return
endif

c Non DDM JOB after restart

if(inc.eq.0.and.iprestt_data(23).ne.0.and.nprocd.eq.0) then
call readresfile
endif

c DDM JOB after restart/Remeshing

if(inc.eq.0.and.iprestt_data(23).ne.0.and.nprocd.gt.0) then
c -- File will be read by each domain, one at a time
do i=1,nprocd
if(i.eq.iprcnm) then
call readresfile
idone=i
endif
call domflag(idone,dummy,1,0,1,0)
enddo !i=1,nprocd
endif
return
end

c -- Function to read user variables

subroutine readresfile
use uservariables
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
#include "op1.cmn"
#include "sfkin.cmn"

integer :: iunit = 18
integer :: kin_unit = 57
character*(MAX_CARD_DIM_2) filen,filen2

call get_kin2res_filename(filen, filen2)

call openfile(iunit,filen,1,0)
call sfkin2res (iunit,0,0,0,0,1,1) !read
integer_array(1:number_integer)=i2res(1:number_integer)
real_array(1:number_real)=x2res(1:number_real)
call closefile(iunit,filen,0,2)
return

67
2023.2 User subroutines Restart with global variables

end

c -- Function to save user variables

subroutine saveuservariables
use uservariables
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
#include "op1.cmn"
#include "sfkin.cmn"
#include "concom.cmn"
#include "cdominfo.cmn"
#include "marc_prestt.cmn"

character*(MAX_CARD_DIM_2) filen,filen2
integer :: iunit = 18
integer :: kin_unit = 57
c Saved by only one domain
if (iparent.eq.0) then
call get_kin2res_filename(filen, filen2)
i2res(1:number_integer)=integer_array(1:number_integer)
x2res(1:number_real)=real_array(1:number_real)
call openfile(iunit,filen,3,0)
call sfkin2res (iunit,kin_unit,
+ number_integer,number_real,0,1,3) !append
call closefile(iunit,filen,0,2)
endif
return
end

c -- Function to get kinematic file name

subroutine get_kin2res_filename(filename1, filename2)
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
#include "op1.cmn"
#include "jname.cmn"
#include "cdominfo.cmn"

integer i1,last_char
character*(MAX_CARD_DIM_2) filename1,filename2

if(nprocd.gt.0)then
i1=last_char(jidnam_base)
filename1=jidnam_base(:i1)//'_kin.res'
filename2=jidnam_base(:i1)//'_kin.res0'
else
i1=last_char(jidnam)
filename1=jidnam(:i1)//'_kin.res'
filename2=jidnam(:i1)//'_kin.res0'
endif
return

68
2023.2 User subroutines Electrical contact coefficient

end

1.5.7. Electrical contact coefficient

With the user subroutine uelcnco(), electrical contact coefficient (conductance) at the interface in a joule heating
simulation can be defined. This becomes useful in process, such as Resistance spot welding, where the electrical
contact coefficient can be a complex function, especially when default electrical contact coefficient calculated by the
solver may not capture the observed behaviour. To use this subroutine the input dat file need to be modified manually.

1.5.7.1. Example
A 2D-Spotwelding process from the demo and example is considered. The provided user subroutine need to be linked
to the simulation through the forming control.

Figure 1.31. Link uelcnco() subroutine in Forming control

1.5.7.1.1. Modify input (*.dat) file

After model setup, the job input data needs to be created and manually modified. Currently the FORMING GUI do
not create the proper flag to support uelcnco() subroutine. The necessary modification is done under the CONTACT
TABLE card entry. Create the process input file and then use the "Edit DAT file" option to open the created dat file
in the default text editor.

User should make sure that the following modifications are carried out for all the contact pairs
entry in the contact table, for which user want the solver to calculate electrical contact conductivity
through the linked subroutine. Also the change need to be made for the same contact pairs under
different loadcases.

Figure 1.32. "Edit dat file" option

Search the dat file for the occurrence of the text "contact table" and look for the contact pair entry similar to following:

69
2023.2 User subroutines Electrical contact coefficient

Consider the 2nd and 3rd line from bottom and change the 4th entry from 1 to 101 so that uelcnco() subroutine is
called for this contact pair.

101 is the equation type which can be used in the user subroutine to implement different formula for different contact
interface. A different equation type for another contact pair can be used to implement different calculation formula in
the subroutine for different contact pair. Any number greater than 100 will trigger solver to call uelcnco() subroutine.
The numbers highlighted as xpar(1), xpar(2) etc. are parameters that can also be changed. The meaning of these are
described in the subroutine headers. Entry of the above two input data cards (lines) can be described as:

Figure 1.33. Contact table entry to support uelcnco() subroutine

Following subroutine is being used in the example. In the example, electrical contact coefficient for all the contact
pairs uses the same formula.

subroutine uelcnco(iequ,mibody,ibody1,ibody2,nodnum,xpar,
& valnod,uecc)
c * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
c user subroutine to returns user electrical contact coefficient
(conductance)
c
c Based on Bay and Wanheim's model described in equation(9.26)

70
2023.2 User subroutines Electrical contact coefficient

c Zhang H, Senkara J. Resistance welding: fundamentals and applications.

(2011)
c * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
c
c INPUT:
c iequ : equation number from contact card
c
c mibody(1) : element number where electrical contact coefficient is
calculated (ibody1)
c (2) : side of said element
c (3) : integration point of said side
c
c ibody1 : body id to which point belongs (touched body, body 1)
c ibody2 : body id that is touching the body "ibody1" (touching body,
body 2)
c nodnum : internal node number of touching node
c
c xpar() : User input parameter from contact table. Values are
extracted from table if table dependent
c (1) : parameter 1 from contact table
c (2) : parameter 2 from contact table
c (3) : parameter 3 from contact table
c (4) : parameter 4 from contact table
c
c valnod() : nodal values of touching node from body 2, unless specified
otherwise
c (1-6) == stresses at nodes
c ( 7) == temperature at nodes
c ( 8) == flow stress at nodes
c ( 9) == previous sliding velocity 1 wrt body 2
c ( 10) == previous sliding velocity 2 wrt body 2
c ( 11) == normal stress at nodes wrt body 2
c ( 12) == sliding velocity 1 wrt body 2
c ( 13) == sliding velocity 2/angle in 2D wrt body 2
c ( 14) == velocity 1 rigid die global wrt body 2
c ( 15) == velocity 2 rigid die global wrt body 2
c ( 16) == velocity 3 rigid die global wrt body 2
c ( 17) == temperature of touched body wrt body 2
c ( 18) == voltage at nodes
c ( 19) == voltage of touched body wrt body 2
c ( 20) == x-component of inward normal wrt body 2
c ( 21) == y-component of inward normal wrt body 2
c ( 22) == z-component of inward normal wrt body 2
c ( 23) == area
c ( 24) == concentration/pressure at nodes
c ( 25) == concentration/pressure of touched body wrt body 2
c ( 26) == temperatures at the nodes begining of increment
c ( 27) == voltage at nodes begining of increment
c ( 28) == concentration/pressure at nodes begining of increment
c ( 29) == temperature of touched body begining of increment wrt body 2
c ( 30) == voltage of touched body begining of increment wrt body 2
c ( 31) == concentration/pressure of touched body begining of inc. wrt body
2
c ( 32) == tangential stress at nodes along slip direction.
c ( 41) == peak temperature
c ( 42) == glue factor
c ( 51) == chemical content at start of increment (at node)
c ( 52) == chemical content change (at node)

71
2023.2 User subroutines FORTRAN- C/C++ user subroutine
interface

c ( 53) == chemical content at start of increment (touched body)

c ( 54) == chemical content change (touched body)
c
c * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
c OUTPUT:
c uecc : user calculated electrical contact coefficient
c * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
c
c
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
c ** Start of generated type statements **
integer iequ,mibody,ibody1,ibody2,nodnum
real*8 xpar,valnod,uecc
c ** End of generated type statements **
dimension mibody(*),xpar(*),valnod(*)
c
real*8 ro1,ro2,ro3,f_stress,sig_n

c xpar(1) = Resistance of body ibody1

c xpar(2) = Resistance of body ibody2
c xpar(3) = Resistance due to contaminant(glue, coating etc.)
c xpar(4) = Thickness of the contact interface

c Calculating resistivity
ro1=xpar(1)*xpar(4)
ro2=xpar(2)*xpar(4)
ro3=xpar(3)*xpar(4)

f_stress=valnod(8) !flow stress

sig_n=abs(valnod(11)) !contact normal stress

if(sig_n.ne.0.0d0) then
c contact resistivity, Bay-Wanheim model
uecc=1.5d0*(f_stress/sig_n)*(ro1+ro2)+ro3

if(uecc.gt.0.0d0) uecc=1.0d0/uecc !conductance

endif

return
end

1.6. FORTRAN- C/C++ user subroutine inter-

face
Simufact Forming supports user subroutines written in C/C++ together with subroutines written in FORTRAN. This
allows the easy integration of already existing C/C++ codes to Simufact Forming with minimal modification.

1.6.1. The interface

The interface works in Windows environment in the following way:

• The user creates a *.dll based on the C/C++ code and puts the *.dll-file in the lib directory of the Simufact installation.

72
2023.2 User subroutines The interface

• Through the FORTRAN subroutine the above created *.dll-file is called through a function. In the argument of the
interface, the user transfers the various variables for the calculation.

A Microsoft Visual Studio project file (sfuserclib.vcxproj) will be provided for easy creation of the *.dll-file from the
C/C++ code. However the *.dll-file can be created by any C/C++ compiler. This tutorial is written for 64bit dll. For
32bit dll changes are required due to handling of the byte sizes for integers and real numbers.

For Linux architecture, the user has to create a dynamic library for linux with the extension .so. With the provided
*.cpp and header file (*.h, referenced in the *.cpp-file), make the library as:

icpc -fPIC -c sfuserclib.cpp

icpc -shared -o libsfuserclib.so sfuserclib.o

The name of the library, libsfuserclib.so, in Linux and sfuserclib.dll in windows are important. The user can link mul-
tiple C/C++ files to create one library. This library should either be put in the lib-directory of the Simufact installation
at C:\Program Files\simufact\forming\2023.2\sfMarc\lib\win64i8 (for windows) or its path
should be added to LD_LIBRARY_PATH (in case of Linux) before submitting the job.

export LD_LIBRARY_PATH="/home/....../path_of_the_library":$LD_LIBRARY_PATH

In FORTRAN the interface is built as a function like:

ret=sfcinterface(t,ni,nd,nc,iarray, darray, name)

where

• ret: An integer value returned from the function. In the absence of the dll, this value will be -1. As C/C++ has an
option to return a value, this argument can be used as a flag to check the successful execution of the code.

• t: An integer value as input. This value can be used to distinguish the different instances when the interface is called.
This value is used to determine which code block to execute in the dll. Through this features the same dll can be
used for various simulations and can be developed simultaneously in time with the added features without breaking
any earlier developed feature. The user just needs to use a different t value.

• ni, nd and nc: The variables corresponding to the number of integers, doubles and characters respectively that will
be passed in the interface.

• iarray, darray and name: Arrays containing integers, doubles and characters information respectively. The size
of these arrays are ni, nd and nc respectively.

The programing interface in the C/C++ file will be like

// sf_c++lib.cpp : Definiert die exportierten Funktionen für die DLL-

Anwendung.
// sf_c++lib.cpp : Defines the exported functions for DLL-application.

#include "sfuserclib.h"
#include <stdio.h>

/*----------------------------------------------------------------
This is the generic interface function for marc calling user c routines.
\param[in] type The type of function that has been called.
\param[in] nbInts The size of the integer array.
\param[in] nbDoubles The size of the double array.
\param[in] nbChar The size of the character array.
\param[in] ints The integer array.
\param[in] doubles The the double array.
\param[in] char The the character array.
------------------------------------------------------------------*/

int64_t sf_cInterface(int64_t * type, int64_t * nbInts, int64_t * nbDoubles,

int64_t * nbChars,

73
2023.2 User subroutines The interface

int64_t * ints, double *doubles, char * chars)

{

switch(*type)
{
case 101:
{
//Test case to print all the infromation of arrays passed from
fortran
//The values will be printed in the jobname.log file
int i;
fprintf(stderr, "Calling default user c interface for function
code %lld\n", *type);
fprintf(stderr, " integer (%lld):\n", *nbInts);
for(i=0; i &lt; *nbInts; i++){
fprintf(stderr, "%d %lld\n",i, ints[i]);
}

fprintf(stderr, " double (%lld):\n", *nbDoubles);

for(int i=0; i &lt; *nbDoubles; i++){
fprintf(stderr, "%d %lf\n",i, doubles[i]);
}

fprintf(stderr, " chars (%lld): '", *nbChars);

fprintf(stderr, "%s", chars);

fprintf(stderr, "'\n");
fprintf(stderr, " returning 101\n");
return 101;
}
break;

case 1: //Case 1
{
// Code related to case 1 comes here
// C/C++ code
//
return 1;
}
break;

case 2: //case 2
{
// Code related to case 2 comes here
// C/C++ code
//
return 2;
}
break;

default: // If type is none of the above

// Do nothing and return unique code, -100 in this example
fprintf(stderr, "Out of bound case, return -100\n");
return -100;
break;

}
}

74
2023.2 User subroutines Example

Some silent features of the interface are:

• The variable values changed in the C/C++ function will also automatically be realized in the FORTRAN side. This
can be used to pass the complete Marc array through the interface. However the user should be careful not to change
the Marc specific array through the interface, because that may cause an error.

• The array length can be changed in the C/C++ side if required. However care should be taken with the indexes as
in FORTRAN array index starts with 1 whereas in C/C++ it starts with 0.

1.6.2. Example
The example subroutine for a user-defined motion is modified to control the maximum force and die velocity through
the C/C++ interface. To run this example, the user should first create and link the dynamic library as explained earlier.

The user subroutine in FORTRAN calls the C/C++ interface with the case flag 1 or 2:

c--------------------------------------------------------------
c This routine demonstrate the control of motion of the die
c during the simulation. Control criteria can be entered
c through GUI which is then used in the subroutine
c
c The motion example routine is modified to get the maximum force
c and velocity from c-dll.
c--------------------------------------------------------------

subroutine ubginc(inc,incsub)

c--------------------------------------------------------------
c Subroutine to get value from c dll
c--------------------------------------------------------------

#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
#include "sf_userdata.cmn" !User variable through gui
#include "machin.cmn" !iprtl
#include "cdominfo.cmn" !ddminfo nprocd

integer inc, incsub, icase, ni,nd,nc,ret,sfcinterface,iarray(1)

integer idummy
real*8 darray(2)
character*(10) cname

real*8 xforce,xvelocity
common/force_vel/xforce,xvelocity
if(inc.eq.0) then
c -- Initial initializing to low value
xforce=0.0
xvelocity=-1.0e19
c -- Do for only one domain
if((nprocd.eq.0).or.(iparent.eq.0)) then
cname='Motion'
ni=1 !no. of integer
nd=2 !no. of double
nc=6 !no. of character
iarray= (/ 0 /)
darray= (/ 0.0,0.0 /)

75
2023.2 User subroutines Example

if(niuser.lt.1) then
write(iprtl,'(/,/7x,a)')
+ '*** Error ***: Missing case variable in gui.'
write(iprtl,'(7x,a,/,7x,a,/)')
+ 'Provide at least 1 integer user data as a case value',
+ 'Simulation halted'
call quit(13)
endif
c -- Get the case id from GUI
icase=iuserdata(1)
c -- Call the c routine with the case flag. Case flag will be used
c -- to seperate the different case call
ret=sfcinterface(icase,ni,nd,nc,iarray, darray, cname)

if(ret.eq.1 .or. ret.eq.2) then

xforce=darray(1) ! Force
xvelocity=darray(2) ! Velocity
write(iprtl,'(a,f8.2,a,f6.2,a)') 'Limiting force = ',
+ xforce/1000,' kN, velocity = ',xvelocity*1000,'mm/sec used.'
endif

if(ret.eq.-1.or.ret.eq.-100) then
if(ret.eq.-1) then
write(iprtl,'(a)') 'sfuserclib.dll does not exist'
elseif(ret.eq.-100) then
write(iprtl,'(a)') 'Unknown test case'
endif

write(iprtl,*)
+ 'Default: Max force=5000KN, velocity = 10mm/s assumed'
xforce=5000*1000
xvelocity=-15*0.001
endif

endif
endif

if(nprocd.gt.0) then
call domflag(idummy,xforce,0,1,0,1)
call domflag(idummy,xvelocity,0,1,0,1)
endif
return
end

c--------------------------------------------------------------
c Subroutine to define die velocity
c--------------------------------------------------------------
subroutine motion(x,f,v,time,dtime,nsurf,inc)

#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif

#include "spaceco.cmn" !contact body information, cbname, numdie_cn

#include "machin.cmn" !iprtl

76
2023.2 User subroutines Example

#include "prepro.cmn" !kou

#include "sf_userdata.cmn" !User variable through gui
#include "cdominfo.cmn" !ddminfo nprocd
#include "dimen.cmn" !ncrd

c ** Start of generated type statements **

integer inc, nsurf
real*8 dtime, f, time, v, x
c ** End of generated type statements **
dimension x(*),v(*),f(*)

character*24 cbody,cbody0,fname
integer i,imvbody,i2switch,ilen1,ifunit
real*8 dummy1,zforce

real*8 xforce,xvelocity
common/force_vel/xforce,xvelocity

logical there

fname="motion.txt"

i2switch=0 !flag for die movement direction: Original(=0), or

reverse(=1)
imvbody=0 ! Body id of the die to monitor the movement

c -- The if command is executed once. In case of ddm job only parent

domain
c -- execute it.
if((nprocd.eq.0).or.(iparent.eq.0)) then
inquire( file=trim(fname), exist=there ) !check if file is
there.
if(.not. there) then
if(nsurf.eq.1) then
i2switch=0
c -- Die name defination of cbody is important.
c -- should be same as defined in gui
cbody='Top-die'
ilen1=index(cbody,' ') - 1
do i=1,numdie_cn
cbody0=cbname(i)
if(cbody0(1:ilen1).eq.cbody(1:ilen1))then
imvbody=i !extracts the body id and saces it
exit !skip remaining do loop
endif
enddo

ifunit=3000
call openfile(ifunit,trim(fname),1,2)
write(ifunit,*) imvbody,i2switch !write the file
call closefile(ifunit,trim(fname),0,2)
endif

else
ifunit=3001
call openfile(ifunit,trim(fname),1,2)
read(ifunit,*) imvbody,i2switch !read from file,
call closefile(ifunit,trim(fname),0,2)

77
2023.2 User subroutines Example

endif

endif

c -- In case of ddm job, share the read value of imvbody & i2switch
c -- Since these balues will be >0 in parent domain and 0 in other domain
c -- Dom flag has been used to share the maximum value of imvbody &
c -- i2switch in all the domains
if(nprocd.gt.0) then
call domflag(i2switch,dummy1,1,0,1,0)
call domflag(imvbody,dummy1,1,0,1,0)
endif

c -- Skip alll bodies which should not move

if(nsurf.ne.imvbody)goto 999

c -- Get the absolute value of z-force on die depending upon 2d & 3d

c -- simulation
if(ncrd.eq.2) then !2D case
zforce=abs(f(1))
else !3d case
zforce=abs(f(3))
endif

c -- In case of original die direction, if force exceed maximum force

c -- then reverse the die direction
if((i2switch.eq.0).and.(zforce.gt.xforce))then
i2switch=1
c -- Done by only parent domain in case of ddm job.
if((nprocd.eq.0).or.(iparent.eq.0)) then
ifunit=3000
c -- Replace the original file with new interval
call openfile(ifunit,trim(fname),1,2)
write(ifunit,*) imvbody,i2switch
call closefile(ifunit,trim(fname),0,2)
endif

c -- Share the new value of i2switch in all the domain

if(nprocd.gt.0) then
call domflag(i2switch,dummy1,1,0,1,0)
endif
endif

c -- If maximum force reached, assign the new die velocity to the die
c -- Velocity sign is important wrt simulation direction
if(ncrd.eq.2) then !2D case
v(1)=xvelocity
if(i2switch.eq.1) v(1) = -v(1)
else !3D case
v(3)=xvelocity
if(i2switch.eq.1) v(3) = -v(3)
endif

999 return
1001 format(13x,'*** error ',
* 'Number of user defined real variables nruser:',i3,'
* is less then needed: ',i3,/13x,'check input deck!')

78
2023.2 User subroutines Example

1002 format(13x,'max force of ',f12.6,' is reached in body:'a,

* ' force is:',f12.6)
1003 format(13x,'Calculated velocity of Body ',a,' is:',f12.6)

end

The example C/C++ interface used to determine the velocity is:

// sf_c++lib.cpp : Definiert die exportierten Funktionen für die DLL-

Anwendung.
//

#include "sfuserclib.h"
#include &lt;stdio.h&gt;

/**
This is the generic interface function for marc calling user c routines.
\param[in] type The type of function that has been called.
\param[in] nbInts The size of the integer array.
\param[in] nbDoubles The size of the double array.
\param[in] nbChar The size of the character array.
\param[in] ints The integer array.
\param[in] doubles The the double array.
\param[in] char The the character array.
*/

int64_t sf_cInterface(int64_t * type, int64_t * nbInts, int64_t * nbDoubles,

int64_t * nbChars,
int64_t * ints, double *doubles, char * chars)
{
switch(*type)
{
case 101:
{
//Test case to print all the infromation of arrays passed from
fortran
//The values will be printed in the job *.log file
int i;
fprintf(stderr, "Calling default user c interface for function
code %d\n", *type);
fprintf(stderr, " integer (%d):\n", *nbInts);
for(i=0; i &lt; *nbInts; i++){
fprintf(stderr, "%d %d\n",i, ints[i]);
}

fprintf(stderr, " double (%d):\n", *nbDoubles);

for(int i=0; i &lt; *nbDoubles; i++){
fprintf(stderr, "%d %lf\n",i, doubles[i]);

fprintf(stderr, " chars (%d): '", *nbChars);

fprintf(stderr, "%s", chars);

fprintf(stderr, "'\n");
fprintf(stderr, " returning 101\n");
return 101;
}
break;

79
2023.2 User subroutines Example

case 1: //motion example 1

{
//max force of 6000KN and velocity of 20 mm/sec converted to
SI Units

doubles[0]=6000*1000;
doubles[1]=-20*0.001;
fprintf(stderr, "Case 1 with force and velocity %lf %lf
\n",doubles[0],doubles[1]);
return 1;
}
break;

case 2: //motion example 2

{
//max force of 4000KN and velocity of 10 mm/sec
doubles[0]=4000*1000;
doubles[1]=-10*0.001;
fprintf(stderr, "Case 2 with force and velocity %lf %lf
\n",doubles[0],doubles[1]);
return 2;
}
break;

default: // If type is none of the above

// Do nothing and return unique code, -100 in this example
fprintf(stderr, "Out of bound case, return -100\n");
return -100;
break;
}
}

The FORTRAN subroutine and the case ID (1 or 2) to determine the maximum force and the velocity from the dll
have to be specified in the forming control under Advanced > User-defined as shown in the picture below:

Figure 1.34. Set the user subroutine and case ID in the forming control

The subroutine umotion will be used, so it has to be selected in the tab User subroutine:

80
2023.2 User subroutines Example

Figure 1.35. Select the subroutine umotion in the forming control

To create the file sfuserclib.dll open the provided file sfuserclib.vcxproj in Microsoft Visual Studio and make sure
that Release and x64 are selected as shown in the figure below:

Figure 1.36. Create the dll in Microsoft Visual Studio

The dll can now be created as shown below:

Figure 1.37. Create the dll in Microsoft Visual Studio

81
Scientific Tutorial
2023.2

2 User defined kinematic

2023.2 User defined kinematic User defined kinematic v1

Educational Objectives

Uses of user subroutine in Simufact Forming for custom process kinematic

Prerequisites

FORTRAN programing knowledge, general programing knowledge

2.1. User defined kinematic v1

User defined kinematic has been developed to model realistic forming steps like automatic reversing rolling processes.
The goal is to automatize the simulation with minimal user input, restricted to input variables, initial geometries and
plan of the rolling passes, to get the required process report from the simulation.

The idea is to extend the current graphical user interface so that the user does not have to make the input in a text
file or in a program source code. A developed template can be used by the user accurately as an integral part of the
software. The implementation is accomplished by developing the interface through an XML file.

User inputs from the graphical interface will be written to a file (.kin), which is then read by the solver and through
the user-subroutine necessary kinematics are implemented in the simulation. After calculation of the kinematics data,
definition of contacts etc., solver calls a writing subroutine to write another file (.kin.hist). Syntax of this file is identical
to the "History Definition" from the Marc Inputs (see Volume C). Solver reads the history data from the .kin.hist file
instead of the original history data written in the input definition of the .dat file. The user subroutine is implemented
in FORTRAN90.

In addition to the previously input data for the process sequence, results from the previous steps are used to create
input data for the next step during the calibration and idle time between the passes. That is, the input variables for the
next pass are calculated during the simulation and for each load step, a new .kin.hist file is written.

The figure below illustrates the relationship between the kinematics and solver. The relevant custom user kinematics
and modules are highlighted in color.

Simufact Forming 15 introduces a new version of the user defined kinematic (see Section 2.3). Both versions have
different functionalities so that the new version does not necessarily replace the previous one. The following chapter
explains the background and basics of the first user defined kinematic version. A lot of this applies to the new version
as well so the contents are important no matter which kinematic version is used.

With Simufact Forming 15 the data format of the XML file that describes the content and appearance of
the user defined kinematic dialog had to be changed. The new file ending is .udk (user defined kinematic
file). It is still an XML file but with a different file suffix. Please change the file suffix of already existing
.xml files to .udk. Additionally all files used by the kinematics file such as images or infosheets need
to be placed in a corresponding sub folder. Sub folder name should be equal to the base name of the
.udk file <baseName>.udk.

83
2023.2 User defined kinematic Graphical interface

Figure 2.1. User kinematic: implementation

With Simufact Forming 15 the definition of language settings dependent translations for the user defined
graphical user interface (GUI) is possible with the use of a so-called translation file (.tr). Their use will
be described in Section 2.2.

This document will step by step explain the process to create and implement the user kinematic for process automati-
zation by incorporating already implemented functionalities.

2.1.1. Graphical interface

The graphical user interface in Simufact Forming is presented in this section. For more information about various
statements in the XML control file are compiled in the section Section 2.1.2.

To understand the feature, first open the standard model for rolling. This is available in the Simufact Demos:

84
2023.2 User defined kinematic Graphical interface

Figure 2.2. Standard rolling model in Simufact Demos&Examples

In the object catalog, new kinematics can be added via right mouse click Press > User-defined.

Figure 2.3. Adding new kinematic in inventory window.

On selection, the following dialog box is opened. In this dialog all files ending with .udk present in the master library
of Simufact Forming and in the folder Additional Files of the current project are displayed. Therefore you can select
between different kinematics. In this view only default user defined kinematic of Simufact Forming can be selected.

85
2023.2 User defined kinematic Graphical interface

Figure 2.4. Dialog box for kinematic selection

On selection, the dialog of this default user defined kinematic is opened:

Figure 2.5. Default user-kinematics example provided in Simufact Forming

It can be used as a template for generating own kinematics. To understand the structure of this dialog box and its at-
tributes, we will use the corresponding template XML file named user-kinematics.udk which is present in the installa-
tion path of Simufact Forming at C:\Program Files\simufact\forming\2023.2\sfForming\Mas-
terlibrary\Kinematics. Before editing, make a copy of the existing file as a backup. When opening the user
defined press as described above, you can now select between these two kinematics. We will modify the file user-
kinematics.udk which must be opened with a text editor in order to start editing. Any editor can be used (e.g. Notepad).
However it is recommended to use an editor which supports syntax highlighting for XML, which simplifies working
on the file. In the following, this file is edited in order to explain the various possibilities.

86
2023.2 User defined kinematic Graphical interface

Figure 2.6. Copy of the Default user-kinematics example for modification

The file user-kinematics.udk consist of various tags. With the tags data can be classified and ordered. A tag starts with
"< >" and is closed by "< />". The tags generally come in pairs, for example:

<dialog name="User-defined kinematics"> [...] </dialog>

The complete dialog tag is described here with [...] between both the dialog tag. Attributes are generally assigned
to the beginning tag, as in this case the name attribute is assigned. Attributes can be changed to assign other properties.

Close the dialog box in Simufact Forming. Then, from the above example change the value in the of the attribute
name from "User-defined kinematics" to "My Name" as

<dialog name="My Name"> [...] </dialog>

Save the modified user-kinematics.udk and then re-open the kinematics as explained above, Press > User-defined,
so now the following dialog opens:

Figure 2.7. Name attribute in the kinematic GUI

The tags can also be used to customize input and processing of data from the user. The complete structure of the file
(without the input field shown by "[...]") is structured as follows:

<?xml version="1.0" encoding="UTF-8"?>

<dialog name="My Name">
<kinematic name="USERKIN1" />
<image filename="user-kinematics.png" />
<allPass>
[Information related to all passes]
</allPass>
<perPass>
[Information related to individual pass]
</perPass>
</dialog>

87
2023.2 User defined kinematic Graphical interface

First line is the header (header data) and is always required. The dialog tag has already been explained above. The
kinematic tag specifies the name of the file to be created, USERKIN1 in this case, and should not be changed unless
there is a compelling reason, such as use of multiple kinematics in the model. The image tag specifies the image to
be displayed in the dialog box. The image should be in the corresponding folder, which is in the same directory as the
user-kinematics.udk and has the same base name. The allPass tag defines the kinematic information to be used in
each pass (in every Heat), together with the information of the available tools used in the model. The perPass tag
contains the information that changes from pass to pass, for example the positioning of the tools or their velocity.

Inside the allPass tag, now three fields are created. One each for pusher, upper roll and lower roll. Therefore the
existing tags and attributes have to be adjusted accordingly. For each pass, also the roll speed and reduction of thickness
between the roll can be specified. The tag-structure creates the information tree which can be visualized as:

dialog

allPass perPass

upper roll roll speed

lower roll
decrease in roll height

manipulator

Figure 2.8. Tree-structure highlighting assignment

of various informations in user-kinematics.udk

The final .udk file after modification of the individual tags looks like this:

Figure 2.9. Screen-shot of the .udk file in a text editor with syntax highlighting

The above modified .udk file when opened in Simufact Forming gives the following dialog:

88
2023.2 User defined kinematic Tags and attributes

Figure 2.10. Dialog box after modifying user-kinematics.udk file

We just generated an individual template that you can use to enter input data for simulation. The exact meaning and
explanation of the tags and their attributes are presented in the following sections.

2.1.2. Tags and attributes

Following are the range of values and input capabilities of the individual attributes:

1. name

The name that is used to identify the input field. It is shown in the GUI if no corresponding label is defined in
a translation file.

2. showOnly

If this attribute is enabled, the field will be grayed out with specific value. The user cannot change the value.

3. isVisible

By default the isVisible attribute is enabled by default even if it is missing. If it is disabled

(isVisible="false"), the value will be hidden in the dialog but is still written into the .kin data file at its
defined kinPos position.

In case the hidden value is inside of the pass table, the corresponding column will be hidden too.

4. dimensionIndex

Dimension values are entered so that the corresponding dimension of the value is determined.

5. unitIndex

Determines which unit system is to be used. If no unitIndex is defined, the SI system is used by default.

89
2023.2 User defined kinematic Tags and attributes

6. valueType

Determines the type of data used.

7. type

Determines the type of field to be generated.

8. singleStep

If the type is spinbox you can define the stepping by which the value is changed each time you click the
corresponding arrow of the spinbox.

9. values

Pre-set with the starting value for the element when the dialog box is called.

10.min

Define a minimum value this input field should have. If the input is below that value it will be highlighted in the
GUI with an orange background. The field type has to be unit or spinbox.

11.max

Define a maximum value this input field should have. If the input is above that value it will be highlighted in the
GUI with an orange background. The field type has to be unit or spinbox.

12.kinPos

Determines the position in which the read data is to be written in the .kin data file. This is an important variable
that is being used both in the .udk file and the user subroutine.

The syntax presented here should always be maintained within a tag. In addition, each value must be properly assigned.
The values specified must be marked with quotation marks at the beginning and end of the value. This means that if
the attribute dimension index should be assigned for pressure, the syntax should be as follows:

dimensionIndex="13"

Comments inside the code can be given by the following tag "<!--" and "-->":

<!-- Comments can be put inside this tag -->

Values assigned to the name attribute and its importance

string Display name for the Field string Display name for the Field

Values assigned to the dimensionIndex attribute and its importance

-1 Dimensionless 0 Acceleration
1 Angle 2 Electric current
3 Energy 4 Force
5 Length 6 Mass
7 Speed 8 Temperature
9 Time 10 Volume
11 Density 12 Power
13 Pressure 14 Thermal conductivity
15 Thermal expansion 16 Absorption (radiation-) dose
17 Heat capacity 18 Specific heat capacity

90
2023.2 User defined kinematic Tags and attributes

19 Force per length 20 Voltage

21 Surface tension 22 Temperature Difference
23 Heat transfer coefficient 24 Strain rate
25 Latent heat 26 Percent
27 Frequency 28 Angular velocity
29 Thermal energy 30 Moment of inertia
31 Heat flux 32 Thermal flux
33 Torque 34 Grain size
35 Unit for data set 36 Torsional stiffness
37 Hardness

Values assigned to the unitIndex attribute and its importance

-1 Without unit 0 SI-unit
1 SI-mm-unit 2 Imperial unit
3 Anglo-American unit 4 Custom unit

Values assigned to the type attribute and its importance

unit Values with unit and dimension
spinbox Values selection (double or integer)
separator Separator line between entries
comboBox Drop-Down list item with tools
checkBox Yes/No Box
comboList Drop-Down list with the user created list
edit Test input field

Values assigned to the valueType attribute and its importance

Tool Tools provided in combo box
string Text input
int Integer input
double Floating point number input

Values assigned to the values attribute and its importance

moved Moving tool
fixed Fixed tool
true/false 1=true and 0=false will be written in .kin file
List From n1:: value1; n2::value2, when value2 is selected, n2 is written in .kin file
Value Value such as "0" - is loaded at the start of the dialogue in the field

Values assigned to the kinPos attribute and its importance

i,j,k i-th data block, j-th value, k-th tool is written in the .kin file

With the various attributes available, it is apparent that not all of them are useful for a particular case. For example,
it is not sensible to connect a tool with a dimension and unit. Such input error will be caught by the program. For
two reasons it is advised to work consistently. First a nonsense entry can confuse the future user and second correct
assignment can remove compatibility issue of the interface with future versions of Simufact Forming.

91
2023.2 User defined kinematic The .kin Data

After inserting explanatory comments, the example looks like this:

Figure 2.11. Screen-shot of modified .udk file with comments

2.1.3. The .kin Data

After entering the input information through the user defined GUI, Simufact Forming creates a .kin file. This data file
contains all the input information. The position of the information in the file is given through the attribute kinPos
which has already been presented. Users with no prior experience with creating a .dat file and unfamiliar to the concept
of data blocks, are advised to check the handbook from Marc manual, "Volume C", where many important topics
are discussed.

In short, the concept is that there are data blocks that contain a certain number of values. It is not mandatory that all
data blocks are always written. It may be possible that a file is complete in three lines with the first line of block 1,
the second of block 2 and a third of block 5. Aim of this section is to explain the complete structure of the .kin file so
that the user can examine if all input data has been read and transferred by Simufact Forming correctly.

The following tables are provided with four columns. The first column is specified with the position (when using the
fixed-format), the second contains the relative position (the first value, the second value ...), the third column contains
the data type (I = integer, A = characters, E = floating point or real numbers). Last column explains the data.

1st data block

Position fix Position free Type Importance
1-10 1. A USERKIN1 (number at the end can vary)

2nd data block

Position fix Position free Type Importance
1-5 1. I Starting increment of the kinematic, default=0 (not supported)
6-10 2. I Kinematic type, default=1
11-15 3. I Error message flag, default=0
16-20 4. I Number of moving and force controlled tool
21-25 5. I Number of stationary object
26-30 6. I Maximum number of integers from the user
31-35 7. I Maximum number of floating point numbers from the user

92
2023.2 User defined kinematic The .kin Data

The third Data block defines the moving tools, power-driven tools at the end of the list. This block is as often repeated,
as defined in the second block at position 4.

3rd data block

Position fix Position free Type Importance
1-24 1. A Name of the moving tool in .dat file

The fourth Data block defines the stationary tools. This block is as often repeated, as defined in the second block at
position 5.

4th data block

Position fix Position free Type Importance
1-24 1. A Name of the stationary tool in .dat file

5th data block

Position fix Position free Type Importance
1-5 1. I Total number of Heat
6-10 2. I Maximum number of pass in one Heat
11-20 3. E Global floating point variable, default=0.0
21-30 4. E Maximum power of the press, default=0.0 (->inactive)

The Blocks 6 and 7 are repeated for every Heat.

6th data block

Position fix Position free Type Importance
1-5 1. I Heat number
6-10 2. I Maximum number of pass in the Heat, default=1

7th data block

Position fix Position free Type Importance
1-5 1. I Heat number
6-15 2. E 1. floating point variable for the actual Heat, default=0.0
16-25 3. E 2. floating point variable for the actual Heat, default=0.0
26-35 4. E 3. floating point variable for the actual Heat, default=0.0
36-45 5. E 4. floating point variable for the actual Heat, default=0.0
46-55 6. E 5. floating point variable for the actual Heat, default=0.0

The blocks 8 to 11 are repeated for every pass in the Heat.

8th data block

Position fix Position free Type Importance
1-5 1. I Heat number
6-10 2. I Pass number in the Heat
11-15 3. I 1. Integer variable for the actual pass, default=0
16-20 4. I 2. Integer variable for the actual pass, default=0
21-25 5. I 3. Integer variable for the actual pass, default=0
26-30 6. I 4. Integer variable for the actual pass, default=0

93
2023.2 User defined kinematic Example of the .kin data

31-35 7. I 5. Integer variable for the actual pass, default=0

36-40 8. I 6. Integer variable for the actual pass, default=0
41-45 9. I 7. Integer variable for the actual pass, default=0
46-50 10. I 8. Integer variable for the actual pass, default=0
51-55 11. I 9. Integer variable for the actual pass, default=0
56-60 12. I 10. Integer variable for the actual pass, default=0

9th data block

Position fix Position free Type Importance
1-5 1. I Heat number
6-10 2. I Pass Number in the Heat
11-20 3. E 1. floating point variable for the actual pass, default=0.0
21-30 4. E 2. floating point variable for the actual pass, default=0.0
31-40 5. E 3. floating point variable for the actual pass, default=0.0
41-50 6. E 4. floating point variable for the actual pass, default=0.0
51-60 7. E 5. floating point variable for the actual pass, default=0.0
61-70 8. E 6. floating point variable for the actual pass, default=0.0
71-80 9. E 7. floating point variable for the actual pass, default=0.0
81-90 10. E 8. floating point variable for the actual pass, default=0.0
91-100 11. E 9. floating point variable for the actual pass, default=0.0
101-110 12. E 10. floating point variable for the actual pass, default=0.0

The Blocks 10 and 11 are also repeated for each movable tools like second block, position 4.

10th data block

Position fix Position free Type Importance
1-10 1. E x-stroke
11-20 3. E y-stroke
21-30 4. E z-stroke
31-40 5. E x-velocity
41-50 6. E y-velocity
51-60 7. E z-velocity

11th data block

Position fix Position free Type Description
1-10 1. E Rotation direction
11-20 3. E Rotational velocity
21-30 4. E x-direction
31-40 5. E y-direction
41-50 6. E z-direction

2.1.4. Example of the .kin data

The structures of the above mentioned table and lists are illustrated in the following example. We consider a rolling
process consisting of two passes and the required input data is entered in the dialog box interface:

94
2023.2 User defined kinematic Example of the .kin data

Figure 2.12. Dialog box with input values

After you enter the values, the dialog box can be closed. The newly created object is then assigned as usual in the
process tree. (If you are not sure how these steps are done, then it is strongly recommended to read the Basic Manual).
You should see the process tree look like this:

Figure 2.13. Process-tree after adding and renaming the kinematic press
The .kin file is created at the beginning of the initialisation process. To just create the .kin file, briefly start the simulation
and cancel the simulation when program moves to "start analysis" (see bottom status bar). The .kin file is created in
the "_Run_" directory of the simulation. To access it simply press the key combination "Ctrl+E" in Simufact Forming

95
2023.2 User defined kinematic Example of the .kin data

window or right click on the process icon and select Open process folder. An explorer window opens with the "_Run_"
directory. The .kin file is inside this directory. It can be opened with any text editor.

As a reminder, here are values of the kinPos attribute used in the .udk file:

• Upper roll:

kinPos="3;1"

• Lower roll:

kinPos="3;2"

• Pusher:

kinPos="3;3"

• Thickness reduction:

kinPos="9;8"

• Roll speed:

kinPos="10;3;1"

The generated .kin file should look like this: (Bold numbers are the entry from the user through the dialog windows)

USERKIN1
$ Start increment, kinematic type, debug level, number of moving tools,
number of fixed tools, number of interger values, number of real values
0,1,0,3,0,10,10,
$ Moving body 1
Upperroll
$ Moving body 2
Lowerroll
$ Moving body 3
Pusher
$ Number of totals superpasses,max number of passes, initial height of dies,
max force for hydr. press
1,2,0.00000E+0,0.00000E+0,
$ Superpass number, number of passes in superpass
1,2,
$ Superpass number, real v1,real v2,real v3,real v4,real v5
1,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,
$ #Superpass,#pass,# of rolling passes,Outputfrequency for moving
pass,Outputfrequency for rolling pass,Outputfrequency for back rolling
pass,user int define 1,user int define 2,user int define 3,user int define
4,user int define 5,user int define 6
1,1,0,0,0,0,0,0,0,0,0,0,
$ #Superpass,#pass,# of rolling passes,factor for step size for moving pass
mani,factor for step size for rolling pass,factor for step size for rolling
back pass,Max. rotation step,user real define 1,user real define 2,user real
define 3,user real define 4,user real define 5
1,1,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E
+0,1.00000E-2,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,
$ Real parameter for Moving Body 1: translation:
xstroke,ystroke,zstroke,xvel,yvel,zvel
0.00000E+0,0.00000E+0,3.35103E+1,0.00000E+0,0.00000E+0,0.00000E+0,
$ rotation: rotdist,rotvel,xaxis,yaxis,zaxis

96
2023.2 User defined kinematic Example of the .kin data

0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,
$ Real parameter for Moving Body 2: translation:
xstroke,ystroke,zstroke,xvel,yvel,zvel
0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,
$ rotation: rotdist,rotvel,xaxis,yaxis,zaxis
0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,
$ Real parameter for Moving Body 3: translation:
xstroke,ystroke,zstroke,xvel,yvel,zvel
0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,
$ rotation: rotdist,rotvel,xaxis,yaxis,zaxis
0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,
$ #Superpass,#pass,# of rolling passes,Outputfrequency for moving
pass,Outputfrequency for rolling pass,Outputfrequency for back rolling
pass,user int define 1,user int define 2,user int define 3,user int define
4,user int define 5,user int define 6
1,2,0,0,0,0,0,0,0,0,0,0,
$ #Superpass,#pass,# of rolling passes,factor for step size for moving pass
mani,factor for step size for rolling pass,factor for step size for rolling
back pass,Max. rotation step,user real define 1,user real define 2,user real
define 3,user real define 4,user real define 5
1,2,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E
+0,1.29999E-2,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,
$ Real parameter for Moving Body 1: translation:
xstroke,ystroke,zstroke,xvel,yvel,zvel
0.00000E+0,0.00000E+0,3.35103E+1,0.00000E+0,0.00000E+0,0.00000E+0,
$ rotation: rotdist,rotvel,xaxis,yaxis,zaxis
0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,
$ Real parameter for Moving Body 2: translation:
xstroke,ystroke,zstroke,xvel,yvel,zvel
0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,
$ rotation: rotdist,rotvel,xaxis,yaxis,zaxis
0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,
$ Real parameter for Moving Body 3: translation:
xstroke,ystroke,zstroke,xvel,yvel,zvel
0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,
$ rotation: rotdist,rotvel,xaxis,yaxis,zaxis
0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,0.00000E+0,

For better understanding of the .kin file, line-by line of the file will be explained. All the lines starting with a "$" are
comment lines and are not read by the system. Also note, that the lines are not terminated by a line feed or line break.
As can be seen from the above output, the numbering does not appear in all the lines, for example from above the
first six printed lines are counted with only 5 lines.

Also note, that the output value from the .kin file does not exactly match with the input value. There are two reasons
for it.

• Rounding errors: e.g. 1.29999E-2 instead of 1.30000E-2

• Unit conversion: The solver works in unit of meter [m] where as the input was in millimeter [mm]. Therefore 13
mm was converted to 1.30000E-2 m in the .kin data (implicit with [m]).

The following table shows the mapping (values that are assigned by the user are in bold). kinPos attribute are shown
only for the first value of the data block.

Allocation table
Line Data block kinPos Comment
1 1 1;1 do not change

97
2023.2 User defined kinematic Solver and subroutine

2 comment
3 2 2;1 do not change
4 comment
5 3 3;1 moving tool 1
6 comment
7 3 3;2 moving tool 2
8 comment
9 3 3;3 moving tool 3
10 comment
11 5 5;1
12 comment
13 6 6;1 6;2 is for number of pass
14 comment
15 7 7;1 7;2 - 7;6 can be used freely (pass independent values)
16 comment
17 8 8;1 is repeated for every pass in a Heat
18 comment
19 9 9;1 9;8 corresponds to "thickness reduction" in Pass 1
20 comment
21 10 10;1;1 10;3;1 corresponds to "roll speed" in Pass 1
22 comment
23 11 11;1;1
24 comment
25 10 10;1;2
26 comment
27 11 11;1;2
28 comment
29 10 10;1;3
30 comment
31 11 11;1;3
32 From here repeat the blocks 8-11 for every additional pass, i.e. all the entries from line 16

2.1.5. Solver and subroutine

2.1.5.1. Concept of sfukin-subroutine
The previous chapter explained the method to enter the user data and store in the file. The stored input data in the .kin
data file is to be used by the sfukin subroutine. First, we review the structure of the job input file (*.dat).

All information regarding the model are presented in the .dat file. This file can be roughly divided into three sections:

1. Parameter section consisting of general information and about the program version and release details.

2. Model definition: consist information about geometries of all the networks, materials and boundaries conditions.
These two sections will not be discussed as they are not relevant to the subroutine concept. These two sections as
well as history definitions are described in Marc manual "Volume C" in details.

98
2023.2 User defined kinematic Solver and subroutine

3. History definition.

The history information is important for the user. All bodies, boundary conditions etc. information are to be used by
the model for next step simulation. Networks and constrain created in the previous load cases has to be considered.
However, addition of new objects or bodies in the "history definition" is not possible. This offers two main advantages:

1. when the value of load cases is not known, all the load cases can be calculated "on-the-fly".

2. load case data is always written after the calculation steps so that when writing a new load case data the information
from previous increment is considered.

.dat Data

Parameters Section

.kin data

Model Definition sfukin-subroutine

(Tools, Networks, Materials, Tabels)

create

replaces
History Definition .kin.Hist data
(Load cases)

Figure 2.14. Concept of sfukin subroutine to define new load case

2.1.5.2. The _kin.hist data

The _kin.hist data can be read similarly as the .kin data. The corresponding data blocks and values are taken from Marc
Manual "Volume C". Not all the maps (as each keywords to be called has to be defined in the data block ) defined in
"Volume C" are defined and activated in the sfukin-subroutine. The essential maps are presented here

• Contact Table: The contact table serves two purposes. First it determines which body is in contact with each other
and what properties this contact have. Second all the bodies shown in the post-processing must be defined inside
the Contact table.

• Motion Change: The "Motion Change" map is used for motion, change of motion (translational and rotation) and
its initiation. It is important to note that a translation or rotation motion will continues as long as until one changes.
That is, if in one case, the load "Motion Change" map is used for a roller to specify a speed, then the "Motion
Change" map in the following load is to be written (with zero entries) to stop the movement. Else the initiated roller
will run as long until it is terminated.

• Release: The "Release" map can be used to release all the bodies in contact. This map should be used together with
the "Move" map, because a contact removed in an increment will be produced again by the contact search algorithm
in the next increment, which is not desirable.

• Approach: The "Approach" map is used to bring the tools in contact with the workpiece. For this to work correctly,
a velocity as well as a search direction must be provided. Therefore, the "approach" map is always used with the
"Motion Change" map.

99
2023.2 User defined kinematic Solver and subroutine

• Move: The "Move" map allows the rigid body motions (translational and rotational). The advantage over the "Mo-
tion Change" map is the fewer required increments to execute it, resulting in relatively short time. This map is
mainly used for positioning of a tool at a defined position.

With these five maps, many kinematics and load cases can be designed. The decisive factor is the combination of
these maps and their configuration, as each map parameter is to be assigned. Before the parameters and variables of
the subroutine sfukin are discussed in following chapters, first we described the method to implement the subroutine
in Simufact Forming.

2.1.5.3. Inside sfukin.f routine

After knowing the process of making the necessary input and linking of solver or external compiler with the Simufact
Forming, the actual subroutine should be programmed. The subroutine is to be written in FORTRAN. Therefore
some general notes on subroutine programming are given next. In the next chapter, arrays and variables are explained
together with the information exchange between the subroutine and the solver.

2.1.5.3.1. General Information

This compiling and linking of the solver is done with the Intel FORTRAN Compiler. The programming of the sub-
routine should be in a fixed format. To implement a functionality in the code in the most efficient way, the following
points should be considered.

1. Concept buildup

Before the actual implementation of a function, a flowchart of the process should be built. This point is important
especially for large system processes. Time invested at the initial planning stage saves time later in the testing phase.

2. Minimum working example

A simple working subroutine should be created first to eliminate errors that may arise due to complex model
structure. A common source of error, which is already discussed in the creation of .kin file, is inconsistency in the
use of measuring units. Therefore the units expected by the solver in the chapter "arrays and variables" should be
kept in mind.

3. Documentation

It is recommended to insert meaningful and extensive comments in the code. This will ensure that the code is
maintainable and is understood by other users. This is especially important (also for the programmer) if the code
is to be used or/and further developed over a long period of time.

4. Testing

After the successful implementation of the concept for a minimal working example, the subroutine should be tested.
This should be tested for many cases as possible. Possible irregularities must be noted and documented so that at
later stage when one changes the parameters and comes across the errors, then the errors can be correctly diagnosed.

5. Model requirements and limitations

Often it is not possible or it requires lot of efforts to handle a difficult emerging problem. For a practical use of the
created subroutine model requirements and limitations must be defined. For example, before the simulation starts,
location of the work-piece and rotation direction of the axis has to be defined.

2.1.5.3.2. Naming convention of variables

Variables names in the subroutine are follow with a specific naming protocols.

• All variables and arrays starting with "x" are of type Real*8 and should be used for real numbers.

• All variables and arrays starting with "i" are of type Integer and should be used for integer numbers.

100
2023.2 User defined kinematic Solver and subroutine

• All variables and arrays starting with "c" are of type Character and should be used to store text information.

• Local variables and arrays should be clearly marked, for example with "uk" after the first character.

So for example to initialize a counter, first character to be used is "i" (integer) followed by "uk" (local variable) and
subsequently with the name of the process to make the counter recognisable. Compliance with these rules and the
explicit definition of all variables (an implicit type definition is strongly discouraged), errors and conflicts can be
avoided. Another advantage is the quick identification of variables and arrays types.

2.1.5.3.3. Switches, counters and parameters

To program the load cases three elements are required:

• Switches

A switch can be set to "1" (= active) or "0" (= inactive). For example, a Move-load case are set to be written, first the
corresponding switch must be set to "1" (in this case i2write (8) = 1). In the routine that calls the subroutine sfukin,
the set switches are evaluated. Before the call sfukin of the switches are set to "0" so that no cards are written by
default. This not apply to the switch always needed, the load case the title, the contact table and the Motion Change .
Write Should not these cards are written, they must be explicitly disabled.

• Counters

As the sfukin subroutine is called in every increment (also for the zero increment) , counters are important to keep
track of the operations. For example, a counter (iuk_count) can be used for loading condition to be activated only
after the third iteration of the routine. The implementation can be like this: (lines beginning with C are comment
lines):

C Check if the increment has reached

if(iuk_count.eq.3)then
C Write the loading condition
xcontrols(44)=1.0d0
elseif(iuk_count.ne.3)then
C Increase counter by one
iuk_count=iuk_count+1
else
call quit(13)
endif

Several things should be clear in the above example . First the counter iuk_count should be initialized. This occurs
generally in the 0-increment (the variable iflag = 0 => 0-increment). Also in the else statement a function "call quit
(13)" is called. With this the program quit with an error message. Generally, it is useful to consider that a situation
may arise when none of the cases is satisfied. Therefore even for the program in cases where the possibility of error
is minimum, it is recommended that the "else" statement is always used to call a "quit" function.

• Parameter

Parameters are needed to specify the boundary conditions. For example to implement a moving boundary condition,
one has to specify which body has to be moved, with or without symmetry, translation, rotation and center of
rotation. An examples is a workpiece with symmetry rotating counterclockwise around the z-axis as a center of
rotation and translated by 10mm along the X-axis can be implemented in the above example as:

C Check if the increment has reached

if(iuk_count.eq.3)then
C Write the loading condition
xcontrols(44)=1.0d0
C One body to be moved
nbd2mv = 1
C Parametering the move of the loading body

101
2023.2 User defined kinematic Solver and subroutine

C Workpiece body number is "1"

xdefbd_move(1,1)=1.0d0
C Symmetry while moving the body
xdefbd_move(2,1)=1.0d0
C 10mm=0.01m translation in the x-axis
xdefbd_move(3,1)=0.01d0
C Rotation around z-axis
C -90°=>clockwise rotation
C +90°=>anti-clockwise rotation
xdefbd_move(8,1)=90.0d0
elseif(iuk_count.ne.3)then
C Increase counter by one
iuk_count=iuk_count+1
else
call quit(13)
endif

For better code readability, it is important to mention the comments and proper units even in small snippets of the
code.

2.1.5.4. Arrays and variables

2.1.5.4.1. Array: xcontrols
The xcontrols array stores various settings and information for the control functions. The details are explained below

xcontrols(9)
Dimension: 9
Data type: real*8
Index I/O Importance/Comment
1 I minimum element length in x-direction divided by 3 (can be use as control parameter)
2 I minimum element length in y-direction divided by 3 (can be use as control parameter)
3 I minimum element length in z-direction divided by 3 (can be use as control parameter)
4 I maximum Rotation per Increment
5 O maximum time for the next step
6 O "control"-Parameter (see volC, control(Histroy) Datablock 3, value 1)
7 O "control"-Parameter (see volC, control(Histroy) Datablock 4, value 1)
8 O Factor for re-meshing, 1.0d0= maintain element; <1.0d0= refine; >1.0d0= coarsen
9 O Process progress [%]

2.1.5.4.2. Array: icontrols

The icontrols integer array stores configuration for various settings. All the values are initialized 0 by default.

icontrols(5)
Dimension: 5
Data type: integer
Index I/O Importance/Comment
1 O Increment number for the next step
2 O "control"-Parameter (see volC, control(Histroy) Datablock 2, value 4)

102
2023.2 User defined kinematic Solver and subroutine

3 O "approach"-Parameter
4 O Write process progress? 0=No : 1=Yes
5 O "load case"-Parameter (Number of active boundary conditions)
6 O Process report parameter, 0=No action, 1=write current increment in post.rep

2.1.5.4.3. Array: i2write

The i2write switch array is used to enable or disable the loading map in the _kin.hist file. By default, 1, 3 and 4 are
initialized to 1.

i2write(10)
Dimension: 10
Data type: integer
Index I/O Importance/Comment
1 O write "title"? 0=No : 1=Yes
2 O write "load case"? 0=No : 1=Yes
3 O write "contact table"? 0=No : 1=Yes
4 O write "motion change"? 0=No 1=Yes
5 O write "release"? 0=No : 1=Yes
6 O write "approach"? 0=No : 1=Yes
7 O write "gravity"? 0=No : 1=Yes
8 O write "move"? 0=No : 1=Yes
9 O Do Remeshing? 0=No : 1=Yes
10 O write stop criteria? 0=No : 1=Yes

2.1.5.4.4. Array: xbd_minmax

This array contains information about the maximum and minimum co-ordinates of the body i in space. By definition,
the first body is the workpiece, followed by the tools (initially moving tools, then fixed) and finally the symmetries.
The values in this array should not be overwritten. It can be used for the control purpose.

xbd_minmax(6,workpiece+ number. of tools (moving and fixed), Number of symmetries)

Dimension: 6*Total number of bodies
Data type: real*8
Index I/O Importance/Comment
1,i I min. x-coordinate of the body i
2,i I max. x-coordinate of the body i
3,i I min. y-coordinate of the body i
4,i I max. y-coordinate of the body i
5,i I min. z-coordinate of the body i
6,i I max. z-coordinate of the body i

2.1.5.4.5. Array: xbd_minmax_atstart

The xbd_minmax_atstart array contains the same entries as the xbd_minmax array but with the initial positions of all
the bodies instead of the current positions.

103
2023.2 User defined kinematic Solver and subroutine

2.1.5.4.6. Array: xbdcntrls

The xbdcntrls array defines the contact and contact conditions as well as some additional properties related to the
individual bodies.

xcontrols(100,Anzahl d. Werkzeuge(bewegt und starr)+Anzahl d. Symmetrien)

Dimension: 100*Number of bodies without workpiece
Data type: real*8
Index I/O Importance/Comment
1,i I/O free
2,i I/O free
3,i I/O free
4,i I/O free
5,i I/O free
6,i I/O free
7,i I/O free
8,i I/O free
9,i I/O free
10,i I/O free
11,i I/O free
12,i I/O free
13,i I/O free
14,i I/O free
15,i I/O free
16,i I/O free
17,i I/O free
18,i I/O free
19,i I/O free
20,i O body i is force controlled? 0.0d0=No : 1.0d0=Yes
21,i - blocked
22,i - blocked
23,i I/O free
24,i I/O free
25,i I/O free
26,i I/O free
27,i I/O free
28,i O change rotational axis of body i? 0.0d0=No : 1.0d0=Yes (Parameter->
xmvbd_currmove(9-11,i))
29,i O change center of rotation of body i? 0.0d0=No : 1.0d0=Yes (Parameter->
xmvbd_currmove(12-14,i))
30,i O body i in contact with workpiece? 0.0d0=No : 1.0d0=Yes
31,i O contact tolerance of body i
32,i O separation force of body i
33,i O Coefficient of friction of body i

104
2023.2 User defined kinematic Solver and subroutine

34,i O Contact bias factor of body i

35,i - blocked
36,i - blocked
37,i - blocked
38,i - blocked
39,i - blocked
40,i - blocked
41,i - blocked
42,i - blocked
43,i - blocked
44,i - blocked
45,i - blocked
46,i - blocked
47,i - blocked
48,i - blocked
49,i - blocked
50,i - blocked
51,i - blocked
52,i - blocked
53,i - blocked
54,i - blocked
55,i - blocked
56,i - blocked
57,i - blocked
58,i - blocked
59,i - blocked
60,i - blocked
61,i - blocked
62,i - blocked
63,i - blocked
64,i - blocked
65,i - blocked
66,i I/O free
67,i I/O free
68,i I/O free
69,i I/O free
70,i O force direction on body i due to Maximum Force Criterion. 1.0d0=x : 2.0d0=y : 3.0d0=z
71,i I/O free
72,i I/O free
73,i I/O free
74,i I/O free
75,i I/O free

105
2023.2 User defined kinematic Solver and subroutine

76,i I/O free

77,i I/O free
78,i I/O free
79,i I/O free
80,i I/O free
81,i I/O free
82,i I/O free
83,i I/O free
84,i I/O free
85,i I/O free
86,i I/O free
87,i I/O free
88,i I/O free
89,i I/O free
80,i I/O free
91,i I/O free
92,i I/O free
93,i I/O free
94,i I/O free
95,i I/O free
96,i I/O free
97,i I/O free
98,i I/O free
99,i I/O free
100,i I/O free

2.1.5.4.7. Array: xdefbd_move

The xdefbd_move array is used as a parameter for boundary condition criteria. Next the variable nbd2mv, determining
the number of bodies to be moved, should be set. This array is a special array where the body index (index:1,i) with
value 1 may be assigned to workpiece as well as tool. Symmetry information tells to follow the applied symmetry
condition while moving. For example xdefbd_move(1,i)=1 represents the work piece if xdefbd_move(2,i)=1 whereas
it represent the tool at kinPos="3;1" when xdefbd_move(2,i)=0.

xdefbd_move(11, Number of moving bodies (must be set through ndb2mv))

Dimension: 11* Number of moving bodies
Data type: real*8
Index I/O Importance/Comment
1,i O body number of the moving bodies (1=Workpiece, 1,...,n=tools as in .kin data)
2,i O moving symmetry? 0.0d0=No : 1.0d0=Yes
3,i O translation in x-direction of the body considered in (1,i) [m]
4,i O translation in y-direction of the body considered in (1,i) [m]
5,i O translation in z-direction of the body considered in (1,i) [m]
6,i O rotation around the x-axis of the body considered in (1,i) [°] (+=counterclockwise)
7,i O rotation around the y-axis of the body considered in (1,i) [°] (+=counterclockwise)

106
2023.2 User defined kinematic Solver and subroutine

8,i O rotation around the z-axis of the body considered in (1,i) [°] (+=counterclockwise)
9,i O x-coordinate of the center of rotation of the body considered in (1,i)
10,i O y-coordinate of the center of rotation of the body considered in (1,i)
11,i O z-coordinate of the center of rotation of the body considered in (1,i)

2.1.5.4.8. Array: xmvbd_param

The xmvbd_param array contains the datablock 10 and 11 of the .kin file. This array make available the information
provided by the user input. It is important that the content of this array, the entered numerical values, must not be
changed in the subroutine.

xmvbd_param(11, number of moving bodies * number of heat * maximum number of passes per Heat)
Dimension: 11* number of moving bodies * number of heat * maximum number of passes per Heat
Data type: real*8
Index I/O Importance/Comment
1,i*j I kinPos=10,1,i in pass j
2,i*j I kinPos=10,2,i in pass j
3,i*j I kinPos=10,3,i in pass j
4,i*j I kinPos=10,4,i in pass j
5,i*j I kinPos=10,5,i in pass j
6,i*j I kinPos=10,6,i in pass j
7,i*j I kinPos=11,1,i in pass j
8,i*j I kinPos=11,2,i in pass j
9,i*j I kinPos=11,3,i in pass j
10,i*j I kinPos=11,4,i in pass j
11,i*j I kinPos=11,5,i in pass j

2.1.5.4.9. Array: xmvbd_currmove

The xmvbd_currmove array is counterpart to the xmvbd_param array providing parameters for the Motion change
load case. If no rotation axis is specified then the predefined rotation of the model is defined.

xmvbd_currmove(14, number of moving bodies * number of heat * maximum number of passes per Heat)
Dimension: 14* number of moving bodies * number of heat * maximum number of passes per Heat
Data type: real*8
Index I/O Importance/Comment
1,i O x-component of translation [m]
2,i O y-component of translation [m]
3,i O z-component of translation [m]
4,i O x-component of velocity [m/s]
5,i O y-component of velocity [m/s]
6,i O z-component of velocity [m/s]]
7,i O length of rotation [rad]
8,i O rotational velocity [rad/s]
9,i O direction cosine of x-component of the axis of rotation [ ] (-> xbdcntrls (28))
10,i O direction cosine of y-component of the axis of rotation [ ] (-> xbdcntrls (28))

107
2023.2 User defined kinematic Solver and subroutine

11,i O direction cosine of z-component of the axis of rotation [ ] (-> xbdcntrls (28))
12,i O x-coordinate of the center of rotation (->xbdcntrls(29))
13,i O y-coordinate of the center of rotation (->xbdcntrls(29))
14,i O z-coordinate of the center of rotation (->xbdcntrls(29))

2.1.5.4.10. Array: xpass_info

The xpass_info array is the data block number 9 from the .kin file. Make sure that content is not changed in the
subroutine.

xpass_info(10, Anzahl der Hitzen * max. Anzahl Stiche pro Hitze)

Dimension: 10, Anzahl der Hitzen * max. Anzahl Stiche pro Hitze
Data type: real*8
Index I/O Importance/Comment
1,i I kinPos=9;3 in pass i
2,i I kinPos=9;4 in pass i
3,i I kinPos=9;5 in pass i
4,i I kinPos=9;6 in pass i
5,i I kinPos=9;7 in pass i
6,i I kinPos=9;8 in pass i
7,i I kinPos=9;9 in pass i
8,i I kinPos=9;10 in pass i
9,i I kinPos=9;11 in pass i
10,i I kinPos=9;12 in pass i

2.1.5.4.11. Array: ipass_info

The ipass_info array stores the information of the 8th datablock of the .kin data file.

ipass_info(10, number of Heat * maximum number of pass per Heat)

Dimension: 10, number of heat * maximum number of pass per Heat
Data type: integer
Index I/O Importance/Comment
1,i I kinPos=8;3 in pass i
2,i I kinPos=8;4 in pass i
3,i I kinPos=8;5 in pass i
4,i I kinPos=8;6 in pass i
5,i I kinPos=8;7 in pass i
6,i I kinPos=8;8 in pass i
7,i I kinPos=8;9 in pass i
8,i I kinPos=8;10 in pass i
9,i I kinPos=8;11 in pass i
10,i I kinPos=8;12 in pass i

2.1.5.4.12. Array: xmachine_info

The xmachine_info array stores the information of the 7th datablock of the .kin data file.

108
2023.2 User defined kinematic Solver and subroutine

xmachine_info(10, number of Heat)

Dimension: 22, number of Heat
Data type: real*8
Index I/O Importance/Comment
1,i I number of pass in the Heat i
2,i I kinPos=7;2 in the Heat i
3,i I kinPos=7;3 in the Heat i
4,i I kinPos=7;4 in the Heat i
5,i I kinPos=7;5 in the Heat i
6,i I kinPos=7;6 in the Heat i
7,i I blocked
8,i I blocked
9,i I blocked
10,i I blocked
11,i I blocked
12,i I blocked
13,i I blocked
14,i I blocked
15,i I blocked
16,i I blocked
17,i I blocked
18,i I blocked
19,i I blocked
20,i I blocked
21,i I blocked
22,i I blocked

2.1.5.4.13. Array: inkininfo

The inkininfo array stores the information of the 2nd datablock of the .kin data file.

ikininfo(7)
Dimension: 7
Data type: integer
Index I/O Importance/Comment
1 I kinPos=2;1 (starting increment of the Kinematics)
2 I kinPos=2;2 (Kinematics type)
3 I kinPos=2;3 (debug level)
4 I kinPos=2;4 (number of moving tools)
5 I kinPos=2;5 (number of static tools)
6 I kinPos=2;6 (number of user editable integer variables)
7 I kinPos=2;7 (number of user editable real variables)
8 - blocked
9 - blocked

109
2023.2 User defined kinematic Example

10 - blocked
11 - blocked
12 - blocked
13 - blocked

2.1.5.4.14. Other variables

The variables which are not listed in the array but used in the subroutine are listed alphabetically in the following
table. In the "type" column "char" stands for character and "int" for the integer types.

Other variables
Dimension: 1
Name Type I/O Importance/Comment
ctitle char O process name, also displayed in the .sts data by Simufact Forming
iflag int I iflag=0 => 0-increment is active
ikintype int I type of kinematics
inc_start int I restart required? 0=No : 1=Yes
load case char O Name of the load case (only for gravity)
mxip int I number of user editable integer variables per pass
mxp int I number of user editable real variables per pass
mxpass int I maximum number of passes in a Heat
mxsp int I number of user editable real variables per Heat
nmbdy int I total number of the bodies in the model
nmfixbdy int I total number of static tools
nmmvbdy int I total number of moving tools
numsuper int I number of Heat

2.1.6. Example
This chapter describes a sample subroutine consisting of examples from the previous chapter. By exploiting the four
fold symmetry of the model, the process consists of the following steps:

1. Due to symmetry, rolling with only top roll and workpiece will be modelled. Thickness reduction will be half of
that is entered.

2. In the first pass, the pusher will move the work-piece towards the roll from left to right.

3. In the second pass, the pusher will move the work-piece towards the roll achieving the two pass rolling from right
to left.

First open the rolling example from the Simufact Demos (Rolling -> Flat Rolling -> Hot Rolling). Four-fold sym-
metry of the roll can be exploited in the model. Following steps are required to prepare the model:

1. Delete the process 'Pass1-2D', 'Pass2-3D' and delete the stage control.

2. Delete the object 'LowerRoller' and 'Pusher-2' in Process 'Pass1-3D'.

3. Delete the object 'Compressed+x' and the contact table in Process 'Pass1-3D'.

4. Delete the press 'Rotation'.

110
2023.2 User defined kinematic Example

5. Add the modified user defined press from the previous sections (right click on presses -> user-defined).

6. Add the press to the process and assign the object 'UpperRoller' and 'Pusher-1' to it.

7. In the object catalogue right click-->remove unused to clean up the catalogue.

8. Translate 'UpperRoller' by 50mm in positive z direction and center the x-position of the roll around the origin.

9. Translate 'Pusher-1' by 100mm in negative x-direction.

10.Translate 'Workpiece' by 50mm in negative x-direction.

11.To speed up the simulation for demonstrati◘on purpose a coarser mesh can be used, e.g. element size =1.8mm.

12.Insert two symmetry planes. First in the x-y plane by z=z(min) and second in the x-z plane by y=y(max) of the
workpiece.

The modified model with the added symmetries should now look like this:

Figure 2.15. Modified model with symmetry planes

The process tree of the model and the object window is shown in the following image.

111
2023.2 User defined kinematic Example

Figure 2.16. Modified process tree and object windows of the model

Since the lower roll is removed, the corresponding entry is deleted from the XML control file. The kinematic dialog
window now looks like (thickness reductions were also adjusted so that a realistic roller case):

Figure 2.17. Modified kinematic dialog window

112
2023.2 User defined kinematic Example

2.1.6.1. Subroutine: sfukin.f

After modifying the model, only source code for sfukin.f has to be written. The user subroutine can be saved in
the folder "Additional files", which is accessible via the GUI or can be found in the project folders of the project.
Additionally the subroutine has to be reference in the forming control. A sample source code for sfukin.f looks like this:

c-------------------------------------------------------------------
c Subroutine for user-defined kinematic example
c Routine is to simulate two pass rolling
c-------------------------------------------------------------------
subroutine sfukin(xmvbd_currmove,xdefbd_move,xmvbd_param,
1 xbd_minmax,xbd_minmax_atstart,xcontrols,nmbdy,nmmvbdy,nmfixbdy,
2 mxsp,mxp,mxip,mxpass,numsuper,nbd2mv,ctitle,iarcoutput,iterm,
3 iflag,xpass_info,ipass_info,xmachine_info,ikininfo,xbdcntrls,
4 loadcase,icontrols,i2write,inc_start,ikintype)

C******************************************************************
C Don't change anything here - these are compiler options *
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
c
#include "sfkin.cmn"
#include "concom.cmn"
#include "machin.cmn"
#include "spaceco.cmn"
C******************************************************************
character*24 ctitle,loadcase

integer nmbdy,nmmvbdy,nmfixbdy,iarcoutput,iterm,iflag,mxsp,
1 mxp,mxpass,numsuper,ipass_info,ikininfo,mxip,ikintype,
2 nbd2mv,icontrols,i2write,inc_start

integer i

real*8 xmvbd_currmove,xmvbd_param,xbd_minmax,
1 xbd_minmax_atstart,xmachine_info,xpass_info,
2 xcontrols,xbdcntrls,xdefbd_move

dimension xmvbd_currmove(14,nmmvbdy),xdefbd_move(11,*),
1 xmvbd_param(14,nmmvbdy*numsuper*mxpass),
2 xbd_minmax(12,nmbdy),xbd_minmax_atstart(12,nmbdy),
3 xmachine_info(mxsp,numsuper),xpass_info(mxp,numsuper*mxpass),
4 ipass_info(mxp,numsuper*mxpass),xcontrols(*),
5 xbdcntrls(100,(nmmvbdy+nmfixbdy+1)),ikininfo(13),
6 icontrols(*),i2write(*)

intent(in) :: xmvbd_param,xbd_minmax,xbd_minmax_atstart,
1 xmachine_info,mxp,mxpass,numsuper,nmfixbdy,iflag,mxsp,
2 xpass_info,ipass_info,ikininfo,nmbdy,nmmvbdy,mxip

C*******************************************************************
C
C These section contains all variables that can be used within
C the sfukin subroutine
C I = Input (parsed in for controlling tasks)

113
2023.2 User defined kinematic Example

C O = Output (variables that can be manipulated by the user)

C # = number
C
C xcontrols(9)
C (1) I -min. element feed size in x-direction /3
C (2) I -min. element feed size in y-direction /3
C (3) I -min. element feed size in z-direction /3
C (4) I -max. way of rotation per increment
C (5) O -max. time for the next step
C (6) O -"control (History)" (-> volC, Block 3, Position 1)
C (7) O -"control (History)" (-> volC, Block 4, Position 1)
C (8) O -remeshing factor (>0 => coarse, <0 =>fine)
C -> activate remeshing with i2write(9)=1
C (9) O - progress of process in %
C -> activate user progress with icontrols(4)=1
C
C icontrols(6)
C (1) O - number of increments for the next step
C (2) O - "control (History)" (-> volC, Block 2, Position 4)
C (3) O - approach parameter
C -> activate approach with i2write(6)=1
C (4) O - write user progress of the process (=1)
C (5) O - load case parameter (number of activated
C boundary conditions)
C (6) O - result flag (a file for an automated process
C report is written, if =1)
C
C i2write(10)
C 1=true 0=false
C (1) O - write title
C (2) O - write load case
C (3) O - write contact table
C (4) O - write motion change
C (5) O - write release
C (6) O - write approach
C (7) O - write gravity
C (8) O - write move
C (9) O - remesh
C (10) O - write terminate criteria
C
C xbd_minmax(6, # of all body)
C - first body is the workpiece (j,1)
C - then moving bodies
C - then fixed bodies
C - then symmetry
C (1,i) I - min. x-coordinate of body i at the moment
C (2,i) I - max. x-coordinate of body i at the moment
C (3,i) I - min. y-coordinate of body i at the moment
C (4,i) I - max. y-coordinate of body i at the moment
C (5,i) I - min. z-coordinate of body i at the moment
C (6,i) I - min. z-coordinate of body i at the moment
C
C xbd_minmax_atstart(6, # of all body)
C - same dimensions and meaning as xbd_minmax
C - position at the beginning of the simulation
C
C xbdcntrls(100, # of all body without the workpiece)
C - unmentioned indices can't be used

114
2023.2 User defined kinematic Example

C (20,i) O - body i is load controlled (1.0d0=true)

C (28,i) O - change rot.axis for body i (1.0d0=true)
C -> directional cos of new axis is given
C in xmvbd_currmove(9-11,i)
C (29,i) O - change center of rotation for body i (1=true)
C -> coordinates of new point are given in
C xmvbd_currmove(12-14,i)
C (30,i) O - body i has contact to the workpiece (1=true)
C (31,i) O - contact tolerance for body i
C (32,i) O - separation force for body i
C (33,i) O - friction coefficient for body i
C (34,i) O - contact bias for body i
C (70,i) O - direction of force for force abort criterion
C 1.0d0=x, 2.0d0=y, 3.0d0=z
C
C xdefbd_move(11,nbd2mv)
C - if xdefbd_move is used nbd2mv has to be set to the number
C of body's that should be used with the move load case
C (1,i) O - body number of the body to move
C (2,i) O - move body with symmetry (1.0d0=true)
C (3,i) O - translation of body in x-direction [m]
C (4,i) O - translation of body in y-direction [m]
C (5,i) O - translation of body in z-direction [m]
C (6,i) O - rotation of body about x-axis [°]
C (7,i) O - rotation of body about y-axis [°]
C (8,i) O - rotation of body about z-axis [°]
C Note: a positive value for the rotation causes
C a counter clockwise rotation
C (9,i) O - center of rotation x-coordinate
C (10,i) O - center of rotation y-coordinate
C (11,i) O - center of rotation z-coordinate
C
C xmvbd_param(11,# of moving bodies * max. # of passes * # of Heat)
C - contains data out of the .kin-file
C (1,i*j) I - kinPos=10;1;i in pass j
C (2,i*j) I - kinPos=10;2;i in pass j
C (3,i*j) I - kinPos=10;3;i in pass j
C (4,i*j) I - kinPos=10;4;i in pass j
C (5,i*j) I - kinPos=10;5;i in pass j
C (6,i*j) I - kinPos=10;6;i in pass j
C (7,i*j) I - kinPos=11;1;i in pass j
C (8,i*j) I - kinPos=11;2;i in pass j
C (9,i*j) I - kinPos=11;3;i in pass j
C (10,i*j) I - kinPos=11;4;i in pass j
C (11,i*j) I - kinPos=11;5;i in pass j
C
C xmvbd_currmove(14,# of moving bodies*max. # of passes * # of Heat)
C - parameter for the motion change load case can be given here
C - all values for body i, where the order of bodies
C is given in the .kin-file
C Note: body 1 is not the workpiece in this case
C (1,i) O - x-component of translation of body i [m]
C (2,i) O - y-component of translation of body i [m]
C (3,i) O - z-component of translation of body i [m]
C (4,i) O - x-component of the velocity of body i [m/s]
C (5,i) O - y-component of the velocity of body i [m/s]
C (6,i) O - z-component of the velocity of body i [m/s]
C (7,i) O - length of the rotational path [rad]

115
2023.2 User defined kinematic Example

C (8,i) O - rotational speed [rad/s]

C (9,i) O - directional cos of x-component of rotaxis []
C (10,i) O - directional cos of y-component of rotaxis []
C (11,i) O - directional cos of z-component of rotaxis []
C Note: activate xbdcntrls(28)
C (12,i) O - x-coordinate center of rotation
C (13,i) O - y-coordinate center of rotation
C (14,i) O - z-coordinate center of rotation
C Note: activate xbdcntrls(29)
C
C xpass_info(10, # of Heat * max. # of passes)
C - contains data out of the .kin-file
C (1,i) I - kinPos=9;3 in pass i
C (2,i) I - kinPos=9;4 in pass i
C (3,i) I - kinPos=9;5 in pass i
C (4,i) I - kinPos=9;6 in pass i
C (5,i) I - kinPos=9;7 in pass i
C (6,i) I - kinPos=9;8 in pass i
C (7,i) I - kinPos=9;9 in pass i
C (8,i) I - kinPos=9;10 in pass i
C (9,i) I - kinPos=9;11 in pass i
C (10,i) I - kinPos=9;12 in pass i
C
C ipass_info(10, # o Heat * max. # of passes)
C - contains data out of the .kin-file
C (1,i) I - kinPos=8;3 in pass i
C (2,i) I - kinPos=8;4 in pass i
C (3,i) I - kinPos=8;5 in pass i
C (4,i) I - kinPos=8;6 in pass i
C (5,i) I - kinPos=8;7 in pass i
C (6,i) I - kinPos=8;8 in pass i
C (7,i) I - kinPos=8;9 in pass i
C (8,i) I - kinPos=8;10 in pass i
C (9,i) I - kinPos=8;11 in pass i
C (10,i) I - kinPos=8;12 in pass i
C
C xmachine_info(10, # of Heat)
C - contains data out of the .kin-file
C - all unmentioned indices can't be used
C (1,i) I - number of passes in Heat i
C (2,i) I - kinPos=7;2 in Heat i
C (3,i) I - kinPos=7;3 in Heat i
C (4,i) I - kinPos=7;4 in Heat i
C (5,i) I - kinPos=7;5 in Heat i
C (6,i) I - kinPos=7;6 in Heat i
C
C ikininfo(7)
C - contains data out of the .kin-file
C - all unmentioned indices can't be used
C (1) I - kinPos=2;1 (starting increment of kinematic)
C (2) I - kinPos=2;2 (kinematic type)
C (3) I - kinPos=2;3 (debug level)
C (4) I - kinPos=2;4 (number of moving tools)
C (5) I - kinPos=2;5 (number of fixed tools)
C (6) I - kinPos=2;6 (number of user defined integer)
C (7) I - kinPos=2;7 (number of user defined real)
C
C other variables

116
2023.2 User defined kinematic Example

C NAME TYPE I/O COMMENT

C ctitle char O process name
C iflag int I iflag=0 => 0-increment is running
C inkintype int I type of kinematic
C inc_start int I if 1 => restart is done
C loadcase char O name of load case (gravity)
C mxip int I max. # of user defined integer per pass
C mxp int I max. # of user defined real per pass
C mxpass int I max. # of passes in a Heat
C mxsp int I max. # of user defined real per Heat
C nmbdy int I # of all bodies in the model
C nmfixbdy int I # of all fixed bodies in the model
C nmmvbdy int I # of all moving bodies in the model
C numsuper int I # of Heats
c iterm: = 1 Flag to indicate the end of the whole job
C
C**********************************************************************
C User defined local variables
integer iuk_operation
c store user defined local variables in local common block
common /sfukin_1/iuk_operation
C**********************************************************************

if(ipass.ne.1) goto 999 ! skip non mechanical pass

if(iflag.eq.0) then
if(inc_start.eq.0) then
C Initializing the counter
iuk_operation=0
endif

if(inc_start.gt.0) then
C reading stored local variable
iuk_operation=i2res(1)
endif

goto 9999
endif

if(iuk_operation.eq.0)then
ctitle='Pass1: Roll positioning'

C Calculating the roll alignment in z direction

C Thickness of workpiece is to be reduced by half
C of the input value (3.0/2==1.5mm), xpass_info(6,1),
C provided in the kinematic-dialog window for the first pass
xdefbd_move(5,1)=xbd_minmax(6,1)-
& 0.5d0*xpass_info(6,1)-xbd_minmax(5,2)

C body number as in the .kin data

xdefbd_move(1,1)=1.0d0
C move the roll without symmetry
xdefbd_move(2,1)=0.0d0
C Only one body has to be moved in this step
nbd2mv = 1
C write the contact table for all the bodies
do i=1,(nmmvbdy+nmfixbdy)
xbdcntrls(30,i) = 1.0d0

117
2023.2 User defined kinematic Example

enddo

C Use small time step for the movement

xcontrols(5)=0.01d0
C rigid body motion by an increment
icontrols(1)=1
C write title, contact table, motion change and movement
i2write(1)=1
i2write(3)=1
i2write(4)=1
i2write(8)=1
C counter increment
iuk_operation=iuk_operation+1

elseif(iuk_operation.eq.1)then
ctitle='Pass1: Workpiece positioning'

C The workpiece has to be moved

xdefbd_move(1,1)=1.0d0
C With Symmetry
xdefbd_move(2,1)=1.0d0
C Positioning the workpiece in x-direction
C with tolerance of 0.0054d0 m
xdefbd_move(3,1)=-xbd_minmax(2,1)
& -sqrt(0.18d0*xpass_info(6,1))+0.0054d0

C One body(work piece) is to be moved in this step

nbd2mv = 1
C writing contact table for all the bodies
do i=1,(nmmvbdy+nmfixbdy)
xbdcntrls(30,i) = 1.0d0
enddo

C small time step for the movement

xcontrols(5)=0.01d0
C rigid body motion by the increment
icontrols(1)=1
C write title, contact table, motion change and movement
i2write(1)=1
i2write(3)=1
i2write(4)=1
i2write(8)=1
C counter increase
iuk_operation=iuk_operation+1

elseif(iuk_operation.eq.2)then
ctitle='Pass1: Pusher Contact'

C Velocity for the searching the contact with the workpiece

C xmvbd_currmove(4,2)==> body with index 2(=pusher) is to
C move with velocity 0.01d0 m/sec in x-direction(index 4)
xmvbd_currmove(4,2)=0.01d0
C writing contact table for all the bodies
do i=1,(nmmvbdy+nmfixbdy)
xbdcntrls(30,i) = 1.0d0
enddo

C Activating the required load control

118
2023.2 User defined kinematic Example

xcontrols(5)=0.01d0
icontrols(1)=1
i2write(1)=1
i2write(3)=1
i2write(4)=1
i2write(6)=1
icontrols(3)=4
C counter increase
iuk_operation=iuk_operation+1

elseif(iuk_operation.eq.3)then
C Actual rolling process takes place in this loop,
C Pusher is pushing the workpiece toward roll for rolling
ctitle='Pass1: Rolling process'
C Translate pusher by 0.01d0 in x-->direction
xmvbd_currmove(1,2)=0.01d0
C Pusher is to be moved with the velocity of
C Convert circumferential speed to linear speed (1/8 factor)
xmvbd_currmove(4,2)=0.0125*3.14d0*0.18d0*xmvbd_param(3,1)
C Roll will rotate by 144 degree (converted to radian)
xmvbd_currmove(7,1)=0.8*3.14d0
C Rotational speed of the roll for first pass,
C taken from the kin file. The value was provided
C by the user in kinematic-dialog window.
xmvbd_currmove(8,1)=xmvbd_param(3,1)
C writing contact table for all the bodies
do i=1,(nmmvbdy+nmfixbdy)
xbdcntrls(30,i) = 1.0d0
enddo

C Activating Required loading maps

i2write(1)=1
i2write(3)=1
i2write(4)=1
C counter increase
iuk_operation=iuk_operation+1

elseif(iuk_operation.eq.4)then
ctitle='Pass2:Preparing'
C In this step roll will be further lowered for second pass
C First body to be move slightly is the workpiece,
C to avoid crushing by the roll during lowering
xdefbd_move(1,1)=1.0d0
C With Symmetry
xdefbd_move(2,1)=1.0d0
C workpiece is moved by 19 mm in x-direction
xdefbd_move(3,1)= 0.019d0

C Second body to be moved is the roll (lowering the roll)

xdefbd_move(1,2)=1.0d0
C With out Symmetry
xdefbd_move(2,2)=0.0d0
C Calculate Roll delivery, in z direction, wrt. first pass
xdefbd_move(5,2)=-0.5d0*xpass_info(6,2)

C Third, the pusher is moved

C Third body is to be moved is pusher with body id=2
xdefbd_move(1,3)=2.0d0

119
2023.2 User defined kinematic Example

C Positioning pusher in x-direction, other side by 250 mm

xdefbd_move(3,3)= 0.25d0

C Total three body is to be moved

nbd2mv = 3
C writing contact table for all the bodies
do i=1,(nmmvbdy+nmfixbdy)
xbdcntrls(30,i) = 1.0d0
enddo

C small time step for the movement

xcontrols(5)=0.01d0
C rigid body motion by the increment
icontrols(1)=1
C write title, contact tabele, motion change and movement
i2write(1)=1
i2write(3)=1
i2write(4)=1
i2write(8)=1
C counter increase
iuk_operation=iuk_operation+1

elseif(iuk_operation.eq.5)then
ctitle='Pass2: Pusher contact'
C Velocity for the searching the contact with the workpiece
xmvbd_currmove(4,2)=-0.01d0
C writing contact table for all the bodies
do i=1,(nmmvbdy+nmfixbdy)
xbdcntrls(30,i) = 1.0d0
enddo

C Activating the required load control

xcontrols(5)=0.01d0
icontrols(1)=1
icontrols(3)=4
i2write(1)=1
i2write(3)=1
i2write(4)=1
i2write(6)=1
icontrols(3)=4

C counter increase
iuk_operation=iuk_operation+1

elseif(iuk_operation.eq.6)then
ctitle='Pass2: Rolling process'
C Translate pusher by 0.01d0 in x-->direction
xmvbd_currmove(1,2)=0.01d0
C Converting the circumferential speed to linear speed
C (with 1/8 factor) in reverse direction
xmvbd_currmove(4,2)=-0.0125*3.14d0*0.18d0*xmvbd_param(3,3)
C Roll will rotate by 144 degree (converted to radian)
xmvbd_currmove(7,1)=0.8*3.14d0
C Rotational speed of the roll for second pass from kin file
C Value provided by user in kinematic window.( - for reverse)

xmvbd_currmove(8,1)=-xmvbd_param(3,3)
C writing contact table for all the bodies

120
2023.2 User defined kinematic Translation files

do i=1,(nmmvbdy+nmfixbdy)
xbdcntrls(30,i) = 1.0d0
enddo

C Activating Required loading maps

i2write(1)=1
i2write(3)=1
i2write(4)=1
C counter increase
iuk_operation=iuk_operation+1

elseif(iuk_operation.eq.7)then
C Process_finished
iterm = 1
else
continue
endif
c
9999 continue
c storing local variable
iresdim=1
i2res(1)=iuk_operation
999 continue
return
end

In the above code, a local variable, iuk_operation, has been used to perform the process stepwise. By default, values of
the local variables are retained throughout the program execution due to compiler settings. This is valid if the process
is run on a single core process. When the program runs with parallel option, the values of local variables may be lost in
between the process steps when there is remeshing. Therefore, the critical local variables should also be saved through
i2res() variable. This can be achieved by including the sfkin.cmn common blocks.

#include "spaceco.cmn"
#include "sfkin.cmn"

2.2. Translation files

With Simufact Forming 2023.2, it is possible to create and use so-called translation files (.tr) in combination with the
user defined kinematic (.udk). With these files you can achieve two things:

1. Use short and universally understandable names for the entries in your .udk (and therefore .kin) file and display
different, more descriptive names in the GUI.

2. Make the more descriptive names language dependent, depending on the language setting used in the Simufact
Forming GUI.

This feature is supported for both the versions of user defined kinematic, i.e. v1 and v2.

2.2.1. Language dependent content

Simufact Forming natively supports multiple language settings. You can select them via Extras > Options > Settings
> General.

121
2023.2 User defined kinematic Creating translation files

Figure 2.18. Languages in Simufact Forming

Each language is assigned a certain suffix. This suffix is added to the filename of the language dependent content. The
following table shows the languages and their corresponding suffixes:

Language Suffix
Chinese (simplified) <filename>_zh-tw
Chinese (traditional) <filename>_zh-cn
German <filename>_de
English <filename>_en
Japanese <filename>_ja
Polish <filename>_pl
Russian <filename>_ru
Turkish <filename>_tr

The logic for the language dependent content is as follows: The GUI will look for the files with the appropriate suffix
depending on the current Simufact Forming language setting. For example if your Simufact Forming GUI is set to Eng-
lish, it will first check for <yourkinematicname>_en.tr (translation file), <yourkinematicname>_en.pdf (infosheet)
and <imagenames>_en.png (images shown in the GUI; can have different file formats). If no such files are available
it will check for the same files (except the translation file) without the suffix. For example if <imagename>_en.png
does not exist it will check for <imagename>.png. The translation file is an exception to this and is only used with
the appropriate suffix.

2.2.2. Creating translation files

The translation file itself has to be created manually. However after the creation, the contents will be generated by
Simufact Forming. Follow these steps to create a translation file:

1. Write your user defined kinematic (.udk file).

2. Create a new file in the corresponding folder (folder has to have the same name as the .udk file and be in the same
directory). The filename should be <yourkinematicname>_<languageSuffix>.tr (e.g. mykinematic_en.tr).

3. Open Simufact Forming and set the language to be the same as <languageSuffix>.

4. Create a new process.

5. Insert your user defined press.

As soon as the the dialog is opened, Simufact Forming will read the translation file. Entries that are missing from
the translation file will be automatically added. Entries that are in the translation file but no longer used in the .udk
file will be marked as obsolete and sorted to the bottom of the translation file.

122
2023.2 User defined kinematic Creating translation files

6. Open the translation file with any text editor. You will see all names used in the user defined dialog with empty
label and description elements.

When creating a new translation file please consider where the current directory is from where the .udk
file is loaded. It is only loaded from the Masterlibrary directory if you open a new project and insert
the press for the first time. As soon as you confirm the press dialog for the first time it will be saved
in the AdditionalFiles folder in this project. If you open this press again in this project it will be loaded
from there and not from the Masterlibrary. In this case a possible translation file needs to be created in
<yourProject>/AdditionalFiles/<yourkinematicname>/.

As stated in point 6, the translation file contains a label and a description entry for every name in the .udk file.
If a translation file exists for the language that Simufact Forming is currently running in, the GUI of the user defined
press will not show the name of the entry used in the .udk file, but the content of the according label in the .tr file.
The content of the description will be shown in the brief description window of the dialog. This allows to use
short names in the .udk file to keep formulas and the .kin file tidy and at the same time use longer, more descriptive
names in the GUI. Simply edit the label and description contents in the translation file as shown in the two
following figures:

Figure 2.19. User defined kinematic file and translation

file side by side (view .kin file in Section 2.3.4)

Figure 2.20. Translated GUI

123
2023.2 User defined kinematic User defined kinematic v2

2.3. User defined kinematic v2

Simufact Forming has a second type of user defined kinematic. This new type allows for the creation of more complex
user defined GUIs and uses a new and improved .kin file format. It can be used, however, to represent virtually any
forming process (which does not rely on pass tables) in Simufact Forming that is not yet covered by official modules.
It was developed to suit the needs of complex ring rolling processes, especially with regard to the use of rolling curves
as input. A simple ring rolling model will be used to demonstrate the new features.

The following chapter explains all the changes that were made. There will be many references to the first kinematic
version so you should be familiar with the previous chapter as well. To use the new features you have to set your
kinematics (.udk) file to version 2 as follows:

<kinematic version="2" [...] />

To use the previous version, simply omit the version attribute.

2.3.1. Graphical interface

Please open the example model from the Demos&Examples browser:

Figure 2.21. User defined kinematic v2 example model

The following figure shows a sample screen shot of the user defined press of the model which was created with the
user defined kinematic v2:

124
2023.2 User defined kinematic Graphical interface

Figure 2.22. GUI features: user defined kinematic v2

The main changes to the graphics introduced with the user defined kinematic v2 are:

Multiple page navigation

One of the main changes is the possibility to add more pages to the user defined dialog. This allows for the input of huge
numbers of input parameters while keeping everything well sorted. The single pages can be grouped in subdivisions
for a better overview. In the .udk file this is defined with the new page tag and subdivision attribute:

<page name="Ring dimensions and finish diameter" subdivision="Ring geometry">

[...] </page>

This is exclusive to user defined kinematic v2.

You can find the .udk file for this ring rolling example in the AdditionalFiles folder of the mentioned
Demos&Examples model.

Editable brief description and linkable infosheet

From Simufact Forming 15 onwards you can use a translation file (.tr) to create a language dependent GUI. In this
translation file you can define additional descriptions for each field which will show in the brief description when you
hover your mouse pointer over the respective field. It is also possible to create language dependent infosheets which
can be directly linked to the user defined dialog. For more information about this please see Section 2.2.

This feature is available for all kinematic versions.

Dialog wide CSV import and export

It is now possible to import or export every value of the GUI from or to a CSV file. This allows for an easy filling
of big dialogs. The format of the CSV file is as follows:

INPUT <!-- marks all regular input variables -->

name,value <!-- the name is defined in the .udk file-->
name,value
.
.

125
2023.2 User defined kinematic Tags and attributes

.
TABLES <!-- marks table input-->
name
x-value1,x-value2,...
y-value1,y-value2,...
name
x-value1,x-value2,...
y-value1,y-value2,...

The name for each value corresponds to the name defined in the .udk file. The import and export always employs the
SI-m unit system. Values which are calculated by a formula will not be exported and cannot be imported. The dialog
wide import does no error checking, so the user is responsible for providing a valid input CSV file.

This feature is only available for user defined kinematic v2.

Removed perPass table

The perPass table is not compatible with user defined kinematic v2 and had to be removed. If you use that feature you
should stick with the previous kinematic version.

2.3.2. Tags and attributes

In addition to the usable tags and attributes in the previous user defined kinematic version (compare Section 2.1.2), v2
offers the possibility to use formulas to automatically calculate input values and to use tables for tabular input data.

Previous and still valid attributes

1. name

The name that is used to identify the input field. It is shown in the GUI if no corresponding label is defined in a
translation file. For kinematic v2, this name should be unique as solver will use this name to differentiate between
the entered values.

2. showOnly

If this attribute is enabled, the field will be grayed out with specific value. The user cannot change the value.

3. isVisible

By default the isVisible attribute is enabled by default even if it is missing. If it is disabled

(isVisible="false"), the value will be hidden in the dialog but is still written into the .kin data file at its
defined kinPos position.

In case the hidden value is inside of the pass table, the corresponding column will be hidden too.

4. dimensionIndex

Dimension values are entered so that the corresponding dimension of the value is determined.

5. unitIndex

Determines which unit system is to be used. If no unitIndex is defined, the SI system is used by default.

6. valueType

Determines the type of data used.

7. type

Determines the type of field to be generated.

8. singleStep

126
2023.2 User defined kinematic Tags and attributes

If the type is spinbox you can define the stepping by which the value is changed each time you click the
corresponding arrow of the spinbox.

9. values

Pre-set with the starting value for the element when the dialog box is called.

10.min

Define a minimum value this input field should have. If the input is below that value it will be highlighted in the
GUI with an orange background. The field type has to be unit or spinbox.

11.max

Define a maximum value this input field should have. If the input is above that value it will be highlighted in the
GUI with an orange background. The field type has to be unit or spinbox.

12.kinPos

Determines the position in which the read data is to be written in the .kin data file. This is an important variable
that is being used in the .udk file.

New attribute

13.formula

Add this attribute to entries with the type of unit and enter a formula to automatically calculate the input field
from other values. Use the name of other entries to use them for the calculation. Their type has to be either a unit
or a spinbox. A formula could look like this:

<entry name="OD_blank" dimensionIndex="5" type="unit" kinPos="3;1"/>

<entry name="ID_blank" dimensionIndex="5" type="unit" kinPos="3;2"/>
<entry name="T_blank" dimensionIndex="5" type="unit" formula="(OD_blank-
ID_blank)/2" kinPos="3;3"/>

Set up like this, the GUI automatically calculates the ring thickness T_blank from the input values of the outer
diameter OD_blank and the inner diameter ID_blank (compare Figure 2.22). Note, that the entry names shown
in the GUI are different because a translation file is used.

The formula syntax is very flexible. For example you can use:

• Arithmetical operators: + , - , * , / , ^

• Functions: sin( ), cos( ), asin( ), sinh( ), lg( ), sqrt( ), abs( ), rad( ), ...

• Constants: pi, e

• Conditional expressions: if, else, elseif

There are numerous checks done to ensure that the used syntax is valid. If the GUI recognizes a formula as invalid
(e. g. uses a variable name that does not exist) it will mark the field with a red background to get the users attention.
You cannot confirm a dialog that contains an invalid formula. The results that are calculated by a valid formula are
set to showOnly="true". This means they are greyed out and the user cannot enter anything in the field.

The use of formulas is exclusive to user defined kinematic v2.

The syntax presented here should always be maintained within a tag. In addition, each value must be properly assigned.
The values specified must be marked with quotation marks at the beginning and end of the value. This means that if
the attribute dimension index should be assigned for pressure, the syntax should be as follows:

dimensionIndex="13"

127
2023.2 User defined kinematic Tags and attributes

Comments inside the code can be given by the following tag "<!--" and "-->":

<!-- Comments can be put inside this tag -->

Values assigned to the name attribute and its importance

string Display name for the Field string Display name for the Field

Values assigned to the dimensionIndex attribute and its importance

-1 Dimensionless 0 Acceleration
1 Angle 2 Electric current
3 Energy 4 Force
5 Length 6 Mass
7 Speed 8 Temperature
9 Time 10 Volume
11 Density 12 Power
13 Pressure 14 Thermal conductivity
15 Thermal expansion 16 Absorption (radiation-) dose
17 Heat capacity 18 Specific heat capacity
19 Force per length 20 Voltage
21 Surface tension 22 Temperature Difference
23 Heat transfer coefficient 24 Strain rate
25 Latent heat 26 Percent
27 Frequency 28 Angular velocity
29 Thermal energy 30 Moment of inertia
31 Heat flux 32 Thermal flux
33 Torque 34 Grain size
35 Unit for data set 36 Torsional stiffness
37 Hardness

Values assigned to the unitIndex attribute and its importance

-1 Without unit 0 SI-unit
1 SI-mm-unit 2 Imperial unit
3 Anglo-American unit 4 Custom unit

Values assigned to the type attribute and its importance

unit Values with unit and dimension
spinbox Values selection (double or integer)
separator Separator line between entries
comboBox Drop-Down list item with tools
checkBox Yes/No Box
comboList Drop-Down list with the user created list
edit Test input field
table Enables the input of tables (details see below)

Values assigned to the valueType attribute and its importance

128
2023.2 User defined kinematic Tags and attributes

Tool Tools provided in combo box

string Text input
int Integer input
double Floating point number input

Values assigned to the values attribute and its importance

moved Moving tool
fixed Fixed tool
true/false 1=true and 0=false will be written in .kin file
List From n1:: value1; n2::value2, when value2 is selected, n2 is written in .kin file
Value Value such as "0" - is loaded at the start of the dialogue in the field

Values assigned to the kinPos attribute and its importance

i,j order in which data is written in kin file, not important for solver perspective

New type

table

You can now add tabular input data to your user defined GUI. The table always consists of 2 columns. The columns
are not associated to a dimension index. This means that you can use these tables for any combination of values, as
long as your user subroutine handles them correctly. You can define default value pairs to be in the table with the
values attribute. The respective code in the .udk file might look like this:

<entry name="RC_MD" valueType="double" type="table" values="0;0; 1;3;"

kinPos="3;6"/>

A table entry in the .udk file will show as a table icon in the GUI. A click on the icon opens the default Simufact
Forming table dialog, where you can manually create a table, view, edit, and im- or export it:

Figure 2.23. Using tabular input data with the user defined GUI
The content of these tables will be written to the .dat file. The .kin file will contain the corresponding table ID.

The type="table" is only available for user defined kinematic v2.

With the various attributes available, it is apparent that not all of them are useful for a particular case. For example,
it is not sensible to connect a tool with a dimension and unit. Such input error will be caught by the program. For

129
2023.2 User defined kinematic The .kin data

two reasons it is advised to work consistently. First a nonsense entry can confuse the future user and second correct
assignment can remove compatibility issue of the interface with future versions of Simufact Forming.

2.3.3. The .kin data

The new possibilities of the user defined kinematic v2 requires an all new, flexible .kin file format. The new format is
easier to understand, read, and is more clearly structured. This chapter will explain the format in detail.

1st data block

Position fix Position free Type Importance
1-10 1. A Enter the word: USERKIN
11-15 2. I Enter version number (2: new kinematic, 0: previous kinematic)
16-20 3. I Enter increment number to terminate analysis (debugging)

2nd data block

Position fix Position free Type Importance
1-5 1. I Start Increment
6-10 2. I Kinematic type
11-15 3. I Debug level
16-20 4. I Total number of variables to read
21-25 5. I Number of integers
26-30 6. I Number of reals
31-35 7. I Number of bodies
35-40 8. I Number of tables

Following block should be repeated for all variables defined in 2nd block, 4th field:

3rd data block

Position fix Position free Type Importance
1-20 1. A Unique name of the variable
21-25 2. I Variable type
1 = Integer
2 = Real
3 = Table ID
4 = Body ID
5 = Combo box index
6 = Check box or option button
7 = Text
26-30 3. I integer variable value (for variable type = 1,5,6)
table ID (for variable type = 3; reserved range: 4500-4999)
body ID (for variable type = 4)
31-40 4. E real variable value (for variable type = 2)
41-80 5. A (Translated) name of variable as shown in GUI
81-120 6. A Body name (for variable type = 4)
Display string (combo box text, for variable type = 5)

130
2023.2 User defined kinematic Example of the .kin data

Text string (for variable type = 7)

There are several requirements that have to be met in your corresponding .udk file. If you make a mistake there,
Simufact Forming will prompt a warning message when opening the press dialog. The requirements are:

• All names have to be unique. However, labels in the translation file (which are shown in the GUI) do not have to be.

• Maximum length of a unique name is 20 characters. See kin file format: 3rd data block, 1st field: A, 1-20.

• Unique names and other strings should not contain '$' (.dat and .kin comment) or ',' (.dat and .kin field separator).

• All entries are stored in kinPos="3;x". You have to use an ascending order starting at "kinPos=3;1",
followed by "kinPos="3;2", "kinPos="3;3", etc. for each page.

• If comboboxes are used to select tool geometries, the values entry must not be fixed. It should be moved or
omitted altogether. No matter if it is actually moved by a press or not.

'$' comments are continuously allowed and will be used to separate the .kin file for easier readability. Every page of
your dialog, as well as every separator line will be automatically added as a comment.

Please see the next chapter for an example .kin file using the new format.

2.3.4. Example of the .kin data

Following is the content of an exemplary .kin file generated from the .udk file that comes with the aforementioned
Demos&Examples model. Please note the different format used between kin file of version 1 and 2. Please compare
the following data with the format table from the previous chapter to see what each individual entry means. As stated
earlier, for tables, the table ID will be written in the .kin file. The actual table will be written to the .dat file under
the corresponding ID.

Also notice that the unique name is the name, defined in the .udk file. The translated name, defined as a label in the
translation file will also be written to the respective entry in the .kin file to allow for a good readability. In comparison
to the previous kinematic version, the order of the entries is more straight forward. Because every entry of every page
starts with kinPos="3;1", they are written to the .kin file in the same order as they are written in the .udk file.
Solver uses the unique name as identifier in this kin format compared to earlier one. Order of variables written in the
kin file is not important for the solver. So, a new variables or page can be added later to the GUI without affecting
the earlier written subroutines.

$ Simufact Forming 2021 (rev 1616162947 x64 (Windows), 2021/03/19)

$
$ Date: 23 Mar 2021 13:29:08
$
USERKIN,2, <- 1st data block
0,0,0,20,0,10,6,2, <- 2nd data block
$
$ Dialog name: 'User defined ring rolling' <- comments
$ Kinematic name: 'USERKIN1'
$
$ Page 1: 'Ring geometry/Ring dimensions and finish..' <- page name
$ Finish dimensions <- separator name
OD_finish,2,0,1.12000E-1,Outer diameter,, <- 3rd data block,
$ repeated for all
$ Page 2: 'Components and strategy/Die association' entries
$ Component Selection
MR,4,6,0.00000E+0,Main roll,MainRoller,
MD,4,7,0.00000E+0,Mandrel,Mandrel,
AX_U,4,3,0.00000E+0,Axial roll upper,Axial-Upper,
AX_B,4,2,0.00000E+0,Axial roll bottom,Axial-Lower,
CR_L,4,4,0.00000E+0,Centering roll 1,CenterRoll1,

131
2023.2 User defined kinematic Solver and subroutine: sfukin2

CR_R,4,5,0.00000E+0,Centering roll 2,CenterRoll2,

$
$ Page 3: 'Components and strategy/Rolling Strategy'
$ Rolling Strategy
AV_MR,2,0,3.14159E-1,Angular velocity of main roll,,
alpha,2,0,3.06305E-1,Half opening angle axial roll ,,
$ Radial Rolling Curve - Ring growth rate vs. ring diameter
RC_MD_RGR,3,4500,0.00000E+0,RGR (y) vs. Diameter (x),KIN_my-kinematic-
v2_RC_MD_RGR,
Z_point,2,0,0.00000E+0,Diameter measuring at z =,,
$ Axial Rolling Curve - Ring thickness / Ring height
RC_AX_U,3,4501,0.00000E+0,Height (y) vs. Thickness (x),KIN_my-kinematic-
v2_RC_AX_U,
$
$ Page 4: 'Components and strategy/Centering rolls geometry'
$ Geometry for centering roll 1
D_CR_L,7,0,0.00000E+0,Diameter,auto,
R_CR_L,2,0,3.85000E-1,Radius,,
H_CR_L,2,0,3.80000E-1,Height,,
L_CR_L,2,0,5.00000E-2,Length,,
$ Geometry for centering roll 2
D_CR_R,7,0,0.00000E+0,Diameter,auto,
R_CR_R,2,0,3.85000E-1,Radius,,
H_CR_R,2,0,3.80000E-1,Height,,
L_CR_R,2,0,5.00000E-2,Length,,

2.3.5. Solver and subroutine: sfukin2

After reading the job input files (.dat, .umt), solver will read the kin file. The first non commented line of kin file
will determine which kinematic control routine to be used for the kinematic simulation. For example, first two line
of a kin file can be:

$ user kinematic
USERKIN,2,

In above example, line beginning with '$' symbol is considered as a comment line and will be ignored by the solver.
In the second line the word "USERKIN,2," informs the solver that the kin file format is for a user-kinematic, version
2. With this information solver will call the kinematic routine sfukin2(). Interface of the sfukin2() is as shown below:

subroutine sfukin2(x,f,v,time,dtime,nsurf,kin_unit,iunit,iflag)
use userkinmodule

c-------------------------------------------------------------
c user subroutine to control user defined kinematic, version 2
c-------------------------------------------------------------
c iunit - unit number for _kin.res (=18)
c kin_unit - unit number for .kin file (=57)
c
c iflag = 0 initialize routine and read data
c called at beginning, increment 0
c also called at beginning after ddm remeshing
c iflag = 1 calculates velocity
c called at every increment, for each cycle and each body (nsurf)
c
c iflag = 2 called at end of increment only once without body information

The sfukin2() routine will be called multiple times during the simulation runtime. These call can be differentiated
by the value of variable iflag.

132
2023.2 User defined kinematic Kin file: reading and saving informa-
tion

• iflag=0: called only once at increment 0 at the beginning of simulation. In this call, rest of the kin file should be read
and saved in the memory. Also, initialization of various runtime flags, memory allocations should be done during
this call. For restart support, this call should be used to read the earlier saved variables.

• iflag=1: called at each cycle during the increment. Within each cycle, the routine is separately called for each body
present in the simulation. So, if the simulation contains 7 bodies, then this routine will be called 7 times with different
body id defined by nsurf. Value of nsurf can be used to know for which body the call is being made. Only during
this call, velocity of the body defined locally will be used by the solver to move the body. Both linear and rotational
velocity can be provided during the call.

• iflag=2: called only once at the end of each increment. This call should be used to output process logs and saving
important process information for next increment or restart.

In case of parallel DDM job, each domain will call the routine separately in parallel. The routine interaction with the
solver is explained by the sketch below

Simufact Forming GUI

Userkin v2 GUI Post

Model setup
( , .TR) processing

Kinematic file Datfile (.dat)

(.kin) Materila(.umt)

Solver

sfukin2 FEM
(.F) Calculation Results

Figure 2.24. sfukin subroutine interaction with the solver

2.3.6. Kin file: reading and saving information

Rest of the kin file should be read during the call with iflag=0 at the beginning of the simulation. The file is already
opened by the solver with unit number kin_unit. The reading of the kin file is done through a call to utility routine

call read_kin_version2(kin_unit,kinerror)

kinerror will be non-zero in case error occurs during reading the kin file. This can be used to stop the simulation and
investigate the cause of error. The information from, 2nd and 3rd data block of kin file is read by above routine and
stored in several variables inside global module name userkinmodule.

For example, the information of second data block from a kin file

0,0,0,30,0,18,6,2,

will be stored by the solver in ival2ndblock() array as

ival2ndblock(1) = 0
ival2ndblock(2) = 0

133
2023.2 User defined kinematic Example

ival2ndblock(3) = 0
ival2ndblock(4) = 30
ival2ndblock(5) = 0
ival2ndblock(6) = 18
ival2ndblock(7) = 6
ival2ndblock(8) = 2

Solver will then use the value of 4th integer values above to allocate the memory for storage of information for 3rd
data block onwards. In above example, solver will allocate various arrays of length 30 and will read the subsequent
kin files. Solver uses six arrays to save the information of all the 3rd block input. For example consider the following
lines from the kin file

D_CR_L,7,0,0.00000E+0,Diameter,auto,
R_CR_L,2,0,3.85000E-1,Radius,,
H_CR_L,2,0,3.80000E-1,Height,,

The above information will be read and stored in various array as shown in figure below.

Array unqname() ivartype() ivarvalue() rvarvalue() aboutvar() aboutvar2()

type character integer integer real character character

Index
i D_CR_L 7 0 0.00E+00 Diameter auto

i+1 R_CR_L 2 0 3.85E-01 Radius

i+2 H_CR_L 2 0 3.80E-01 Height

Figure 2.25. 3rd data block storage variables and type

The content of kin file will be accessible by including the module as

subroutine sfukin2(x,f,v,time,dtime,nsurf,kin_unit,iunit,iflag)
use userkinmodule

Unique variable name assigned in xml file is used to identify the array index information. A utility routine, index_info,
can be used to get the array index information of D_CR_L in above example. With the array information.

idx_height1=index_info('H_CR_L')

The extracted index information, idx_height1, should be saved in a common block so that the search is done only
once and not repeated every time. Also note the search is case insensitive. All the 3rd block informations from kin
file will be available in subroutine as

unqname(idx_height1) will contain "H_CR_L"

ivartype(idx_height1) will contain "2"
ivarvalue(idx_height1) will contain "0"
rvarvalue(idx_height1) will contain "3.8000E-1"
aboutvar(idx_height1) will contain "Height"
aboutvar2(idx_height1) will contain ""

2.3.7. Example
To view the example please open the User kinematics 2: ring rolling model from the Demos&Examples browser
in the scientific category:

134
2023.2 User defined kinematic Example

Figure 2.26. Open the User kinematics 2: ring rolling model

2.3.7.1. Important files for this example

This model contains a RAW process with simplified kinematics. The machine input is done via a user-defined GUI
and the die movements are controlled by a included user subroutine. You will find all the important files (.udk and .f)
in the Additional files folder of the project:

Figure 2.27. Important files in the Additional files folder

The files used to complete the user-defined GUI are in the my-kinematic-v2 folder in Additional Files. Here you
find the images used, as well as the translation files for English and German each. You should take some time to see
if you understand the structure of the .udk file and what is shown in the GUI.

135
2023.2 User defined kinematic Example

Figure 2.28. Comparison between user-defined GUI and respective .udk file

As explained in the previous chapter, it is possible to show translated names (labels) in the GUI instead of the unique
names used in the .udk file (and transferred to the .f subroutine via the .kin file). The respective .tr translation files
which control that are also located in the Additional files folder. Following a comparison between the English .tr file
and the .udk file:

Figure 2.29. Comparison between .tr and .udk file

2.3.7.2. Capabilities and limitations of this ring rolling model

We provide all necessary files for this model which are needed to run the simulation. They can also be edited by
experienced users to change the user-defined GUI or the control algorithms. As real ring rolling control algorithms are
very complex and well guarded secrets by the machine manufacturers, we provide only very basic control algorithms
for this demo model. Before reading through the complete subroutine, we will summarize what this ring rolling control
does and what it does not.

136
2023.2 User defined kinematic Example

• Required input parameters

• Finish dimension - Outer diameter: Acts as the single termination criteria. When this ring diameter is reached,
the simulation stops.

• Die association: Select which die component in the process covers which role of the ring rolling tools.

• Angular velocity of main roll: Basis for the calculation of the axial roll velocities.

• Axial roll angle: Required for the calculation of the axial roll rotational velocities.

• Radial and axial rolling curves: Required table input to control radial and axial reduction as functions of current
ring diameter or ring thickness.

• Centering rolls geometries: Required to calculate the circular path on which the centering rolls will move.

• Implemented die kinematics

• Main roll

• Constant angular velocity

• Mandrel

• Rotating freely

• Radial reduction / Y-velocity based on table input of ring growth velocity

• Axial rolls

• Angular velocity calculated to be slip-free at the outer diameter of the ring

• Axial frame movement /Y-velocity is synchronized with ring growth (ring position between axial rolls is al-
ways constant)

• Synchronized between upper and lower axial roll

• Axial reduction / Z-velocity based on table input

• Centering rolls

• Purely geometrical centering roll positioning assuming an ideal ring geometry

• Rotating freely

You can change the die geometries to adapt this model to your own ring rolling process. If you want to use profiled
geometries you will have to adapt the control routine to account for this change. You can additionally implement your
own control algorithm to make the subroutine more realistic. The subroutine supports DDM and restart functionality.
This has to be considered, when adapting it to your need. Please note that for any changes in the subroutine, a compiler
is needed, to update the solver executable. The compilation of the subroutine has to be specified in the Forming Control.

2.3.7.3. Subroutine: sfukin2.f

The user subroutine used for this model, which is also located in the Additional files folder. The commentary inside
the subroutine should help you understand the code. The routine can be divided into three sections based on the value
of iflag when called.

1. iflag=0

The call will be done once at the beginning at zero increment. In the provided example, following actions are done.

• defining the custom loadcase name

137
2023.2 User defined kinematic Example

• reading rest of the kin file

• extracting the index information of important input variables and saving in local common block for later use

• calculating other fixed parameters like axis of centre roll

• Memory allocation for bodies that should freely rotate based on the contact. Die insert should be used in the model
setting for these bodies. After memory allocation, the control node of translation, jfnold(1,bodyid), is saved in
the kinbcnd for the motion. The velocity of these bodies should then be applied through apply_fixed_bc()
routine.

2. iflag=1

The call is made at every cycle within an increment. Furthermore the call will be made for each body (nsurf) in the
simulation. Analysing the nsurf, the following actions are done.

• If called for workpiece, calculations based on the rolling curves which are time dependent are done and stored
for later use.

• First the outer diameter is calculated by averaging the bounding box dimensions. Additionally the workpiece
center coordinates are stored.

• The positioning of the axial rolls with respect to the workpiece position is stored, to manage the axial rolls y-
movement based on the ring growth. If the ring diameter increases, the axial rolls follow this movement.

• Rotational velocity of axial rolls is calculated

• Ring thickness is calculated based on averaging the two workpiece thicknesses measured at the x=0 plane at
the specified z-value.

• Z-Velocity of upper axial rolls is calculated. The target height for the calculated thickness is extracted from
the table input height-thickness curve. Additionally the current height is calculated as the distance of the axial
rolls. If the height decreases, the z-velocity is set such that the height difference is achieved in this increment.
Note that the z-velocity is always negative, so the axial rolls cannot move back again.

• Calculate mandrel velocity. This is done by evaluating the ring growth velocity curve for the current diameter.
If the ring growth velocity is positive, the target diameter is calculated and the mandrel feed is determined
based on current ring thickness and a rough volume estimation.

• If called for main roll, the angular velocity from input is provided for main roll rotation.

• If called for mandrel, the calculated y-velocity (see call workpiece) is applied . As the mandrel needs free rotation,
the y-velocity needs to be applied through the fixed boundary conditions.

• If called for bottom axial roll, y- and rotational velocity is provided, which were calculated when call was made
for workpiece.

• If called for upper axial roll, the previously calculated y- and rotational velocity as well as the z-velocity is
applied.

• If called for center roll 1 or 2, the roll is moved to the target position considering intersection of two circles
for calculating the center. The velocity is applied through the application of fixed boundary conditions for free
rotation.

3. iflag=2

The call is made once at end of each increment. In the provided example following actions are done in this call:

• ring diameter is calculated to be updated for THS plot. Value of xkinbc(1) and xkinbc(2) from sfkin.cmn is
updated with current and previous ring diameter respectively.
138
2023.2 User defined kinematic Example

• log output is written about ring profile in the out file.

• current ring diameter is checked against the desired ring diameter. If the check is positive, then simulation is
stopped by setting the variable iotnum from iautcr.cmn to next increment number. Simulation will then not
continue in the next increment.

The subroutine used in the example is listed verbatim below. Several Marc specific variables have been used inside.
All the variables that have not been locally defined will be defined in the common blocs. To know more information
about it, please search the variables name against the list of common blocks under Simufact Forming installation
directory at C:\Program Files\simufact\forming\2023.2\sfMarc\common

module sfukin2_module
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
real*8 cr_l_armx, cr_l_army, cr_r_armx, cr_r_army
real*8 tang_vel_mr, term1, yoffset_ini
endmodule sfukin2_module

subroutine sfukin2(x,f,v,time,dtime,nsurf,kin_unit,iunit,iflag)
use userkinmodule
use sfukin2_module

c-------------------------------------------------------------
c user subroutine to control user defined kinematic, version 2
c-------------------------------------------------------------
c iunit - unit number for _kin.res (=18)
c kin_unit - unit number for .kin file (=57)
c
c isfkin=15: User defined kinematic version 2
c
c iflag = 0 initialize routine and read data
c called at beginning, increment 0
c also called at beginning after ddm remeshing
c iflag = 1 calculates velocity
c called at every increment, for each cycle and each body(nsurf)
c
c iflag = 2 called at end of increment only once without body information
c
c
c The below value is only true in case when call is made with iflag=1
c
c input : nsurf - surface(body) number for which data is
requested
c time - the time at which data is requested
c dtime - the current time increment
c x(6) - current die defining coordinates:
c x(1) = 1st coordinate of center of
c rotation
c x(2) = 2nd coordinate of center of
c rotation
c x(3) = 3rd coordinate of center of
c rotation
c Axis for specifying angular velocity:
c x(4) = 1st component of direction cosine
c x(5) = 2nd component of direction cosine

139
2023.2 User defined kinematic Example

c x(6) = 3rd component of direction cosine

c
c f(6) - the current surface load:
c f(1) = 1st component of load
c f(2) = 2nd component of load
c f(3) = 3nd component of load
c f(4) = 1st component of moment
c f(5) = 2nd component of moment
c f(6) = 3rd component of moment
c
c output : v(4) - current surface velocities
c v(1) = 1st component of the velocity at
c the center of rotation
c v(2) = 2nd component of the velocity at
c the center of rotation
c v(3) = 3nd component of the velocity at
c the center of rotation
c v(4) = angular velocity around axis defined
c above with x(4), x(5), and x(6).
c
c

#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
#include "sfkin.cmn"
#include "concom.cmn"
#include "iautcr.cmn"
#include "machin.cmn"
#include "memory.cmn"
#include "prepro.cmn"
#include "feature.cmn"
#include "spaceco.cmn"
#include "creeps.cmn"
c ** Start of generated type statements **
c-------------------------------------------------------
integer nsurf,kin_unit,iunit,iflag
real*8 dtime,f,time, v,x
dimension x(*),v(*),f(*)
c-------------------------------------------------------
c ** End of generated type statements **

integer i,index_info,nintersec,kinerror
integer idx_height1,idx_height2,idx_length1,idx_length2
integer idx_mainroll,idx_mandrel,idx_ax_upper,idx_ax_bottom
integer idx_cr_left,idx_cr_right,idx_md_rgr_tab
integer idx_axu_tab,idx_larm_rad,idx_rarm_rad
integer idx_od,idx_id,idx_wt,idx_ht

real*8 xmin(3),xmax(3),x_mr,y_mr,yy(4)
real*8 r_thick,ring_height,axr_rvel,axr_yvel,axr_zvel,h_up_axr
real*8 r_ht,z_low,z_up,r_ht_curr,radius_mr
real*8 cr_x,cr_y,cr_r,sx1,sy1,sx2,sy2,sin_alpha,yoffset
real*8 roll_x,roll_y,arm_radius,dia_inner,half_angle
real*8 ring_dia,wp_x,wp_y,wp_r
real*8 half,z_measure,yvel_mand,rgv,d1,rt1,rt2,dx,v0

140
2023.2 User defined kinematic Example

character*20 name_ring_rolling
logical calc_ring_ypos

parameter(name_ring_rolling='V2 ring rolling')

parameter(half=0.5d0)

common/my_cmn/idx_mainroll,idx_mandrel,idx_ax_upper,kinerror,
& idx_ax_bottom,idx_cr_left,idx_cr_right,z_measure,
& idx_md_rgr_tab,idx_od,idx_id,idx_wt,idx_ht,
& idx_axu_tab,wp_x,wp_y,wp_r,idx_larm_rad,idx_rarm_rad,
& sin_alpha,axr_rvel,axr_yvel,axr_zvel,yvel_mand

c=======================================================================
if(iflag.eq.0) then
c-----------------------------------------------------------------------
c iflag = 0 Initialization
c-----------------------------------------------------------------------
ifeat(195)=0
call change_lcase_title(trim(name_ring_rolling))
call read_kin_version2(kin_unit,kinerror)

c getting index information from kin input

idx_mainroll=index_info('MR')
idx_mandrel=index_info('md')
idx_ax_upper=index_info('AX_U')
idx_ax_bottom=index_info('AX_B')
idx_cr_left=index_info('CR_L')
idx_cr_right=index_info('CR_R')

idx_md_rgr_tab=index_info('RC_MD_RGR')
idx_axu_tab=index_info('RC_AX_U')

idx_larm_rad=index_info('R_CR_L')
idx_rarm_rad=index_info('R_CR_R')

idx_od=index_info('OD_finish')
idx_id=index_info('ID_finish')
idx_wt=index_info('T_finish')
idx_ht=index_info('H_finish')

c calculation of centre roll's arm hinge point

c Mainroll mid point(x,y) by bounding box
call sfkinbd(ivarvalue(idx_mainroll),xmin,xmax,0)
x_mr=half*(xmin(1) + xmax(1))
y_mr=half*(xmin(2) + xmax(2))

radius_mr=half*(xmax(2)-xmin(2))

idx_height1=index_info('H_CR_L')
idx_height2=index_info('H_CR_R')

idx_length1=index_info('L_CR_L')
idx_length2=index_info('L_CR_R')

c Calculate arm axis of left center roll

cr_l_armx=x_mr + rvarvalue(idx_height1)
cr_l_army=y_mr + rvarvalue(idx_length1)

141
2023.2 User defined kinematic Example

c Calculate arm axis of right center roll

cr_r_armx=x_mr - rvarvalue(idx_height2)
cr_r_army=y_mr + rvarvalue(idx_length2)

c output to jobname.log file once

write(iprtl,'(/)')
write(iprtl,8001)'centre of main roll ',x_mr,y_mr
write(iprtl,8001)'CR left arm axis ',cr_l_armx,cr_l_army
write(iprtl,8001)'CR right arm axis ',cr_r_armx,cr_r_army
write(iprtl,'(/)')

c Memory allocation using internal marc routine for those bodies

c with die insert
mxkinbcnd=3 !1=Mandrell, 2= left CR, 3=right CR
call marcfalloc(mxkinbcnd*3,8,1,pxkinbcvel,
& 'pxkinbcvel ',mem_none,kinerror)

call marcfalloc(mxkinbcnd*3,8,1,pxkinbcdsp,
& 'pxkinbcdsp ',mem_none,kinerror)

call marcfalloc(mxkinbcnd*3,4,1,pkinbcnd,
& 'pkinbcnd ',mem_none,kinerror)

c Saving control nodes for the load controlled bodies

do i=1,3
kinbcnd(i,1)=jfnold(1,ivarvalue(idx_mandrel))
kinbcnd(i,2)=jfnold(1,ivarvalue(idx_cr_left))
kinbcnd(i,3)=jfnold(1,ivarvalue(idx_cr_right))
enddo

c tangential velocity of Main roll

tang_vel_mr=radius_mr*rvarvalue(index_info('AV_MR'))
c save the axial role geometric properties for its velocity calculation
c half angle of conical centre roll
half_angle=rvarvalue(index_info('alpha'))
if(half_angle.eq.0.0d0) then
write(kou,7001)half_angle
write(iprtl,7001)half_angle
kinerror=kinerror+1
endif

c upper axial roll z thickness

call sfkinbd(ivarvalue(idx_ax_upper),xmin,xmax,0)
h_up_axr=xmax(3)-xmin(3)

yoffset_ini=xmin(2)

c saving following value for later

sin_alpha=dsin(half_angle)
term1=h_up_axr*dsin(half_angle)/dsin(2.0d0*half_angle)

c initial leftmost offset between workpiece & upper axial roll

c this offset will be used to move axial roll in -y direction
call sfkinbd(1,xmin,xmax,0)
yoffset_ini=xmin(2)-yoffset_ini

c Initial cr rotational velocity

142
2023.2 User defined kinematic Example

axr_rvel=tang_vel_mr/(term1-sin_alpha*yoffset_ini)

z_measure=rvarvalue(index_info('Z_point'))

if(z_measure.gt.xmax(3) .or. z_measure.lt.xmin(3)) then

write(iprtl,8007) z_measure,xmin(3),xmax(3)
z_measure=half*(xmin(3) + xmax(3))
endif

write(iprtl,8008)z_measure
write(kou,8008)z_measure
call change_lcase_title(trim(name_ring_rolling))
call sfukin2_restart(iflag)

c=======================================================================
elseif(iflag.eq.1) then
c-----------------------------------------------------------------------
c iflag = 1 velocity calculation
c-----------------------------------------------------------------------
c skipping calculation for non mechanical pass and with
c zero time increment

if(dtime.ne.0.0d0 .and. ipass.eq.1 .and.inc.gt.0) then

call change_lcase_title(trim(name_ring_rolling))
kinerror=0
v(1:4)=0.0d0

c if called for workpiece.

if(nsurf.eq.1) then
c various velocities are calculated in this call
c
c calculate ring dia to be used later
call sfkinbd(1,xmin,xmax,0)
ring_dia=half*(xmax(2)+xmax(1)-xmin(2)-xmin(1))

c save workpiece centre coordinate and RADIUS

wp_x=half*(xmax(1) + xmin(1))
wp_y=half*(xmax(2) + xmin(2))
wp_r=0.25d0*(xmax(1)-xmin(1) + xmax(2)-xmin(2))

c calculate y(axr_yvel) and rotational velocity (axr_rvel)

c of axial rolls
c offset between upper axial roll and workpiece
yoffset=xmin(2) !ring left most y value
call sfkinbd(ivarvalue(idx_ax_upper),xmin,xmax,0)
yoffset=yoffset-xmin(2)

axr_yvel=0.0d0
if(yoffset.lt.yoffset_ini) then
axr_yvel=(yoffset-yoffset_ini)/dtime
endif
axr_rvel=tang_vel_mr/(term1-sin_alpha*yoffset)

c calculate z(axr_zvel) y velocity for upper axial roll

c calculate ring inner thickness at z=0.0d0
r_ht=0.0d0
axr_zvel=0.0d0
if(calc_ring_ypos(z_measure,yy,1)) then

143
2023.2 User defined kinematic Example

rt1=yy(2)-yy(1) !left
rt2=yy(4)-yy(3) !right
r_thick=half*(rt1+rt2) !average of left and right

c desired ring height from table

call get_val_table(ivarvalue(idx_axu_tab),r_thick,r_ht)

c calculate ring height based on axial roll positions

call sfkinbd(ivarvalue(idx_ax_bottom),xmin,xmax,0)
z_low=xmax(3)

call sfkinbd(ivarvalue(idx_ax_upper),xmin,xmax,0)
z_up=xmin(3)

r_ht_curr=z_up-z_low !ring height

if(r_ht .lt.r_ht_curr) then !move upper roll in -z direction

axr_zvel=-(r_ht_curr-r_ht)/dtime
else
r_ht=r_ht_curr
endif

c get Ring growth velocity as a function of ring diameter

rgv=0.0d0
call get_val_table(ivarvalue(idx_md_rgr_tab),ring_dia,rgv)
if(rgv.gt.0.0d0) then
d1=ring_dia+rgv*timinc
c rough voulume calculation with cylinderical ring assumption
v0=((ring_dia)**2 -(ring_dia-rt1-rt2)**2)*r_ht_curr/r_ht
dx=rt1+rt2+sqrt(d1*d1-v0)-d1
if(dx.lt.0.0d0) dx=0.0d0
yvel_mand=dx/dtime
endif
endif
endif

c if called for main roll

if(nsurf.eq.ivarvalue(idx_mainroll)) then
v(4)=rvarvalue(index_info('AV_MR')) ! angular velocity of
mainroll
endif

c if called for mandrel

if(nsurf.eq.ivarvalue(idx_mandrel)) then
v(2)=yvel_mand
call apply_fixed_bc(nsurf,1,v(1:3))
endif

c if called for axial roll bottom

if(nsurf.eq.ivarvalue(idx_ax_bottom)) then
v(2)=axr_yvel
v(4)=axr_rvel
endif

c if called for axial roll upper

if(nsurf.eq.ivarvalue(idx_ax_upper)) then
v(2)=axr_yvel
v(3)=axr_zvel

144
2023.2 User defined kinematic Example

v(4)=axr_rvel
endif

c if called for centre roll left

if(nsurf.eq.ivarvalue(idx_cr_left)) then
c curernt position of roll

call sfkinbd(ivarvalue(idx_cr_left),xmin,xmax,0)
c Calculate roll center coordinates
cr_x=half*(xmax(1) + xmin(1)) !x
cr_y=half*(xmax(2) + xmin(2)) !y
cr_r=half*(xmax(1) - xmin(1)) !radius

arm_radius=rvarvalue(idx_larm_rad)

c Get the position of two circle intersection

call twocirc_intersec(0.0d0,wp_y,(wp_r+cr_r),
& cr_l_armx,cr_l_army,arm_radius,
& sx1,sy1,sx2,sy2,nintersec)
if (nintersec.eq.2) then
c find the "left" intersection point
if (sx1.gt.sx2) then
roll_x = sx1
roll_y = sy1
else
roll_x = sx2
roll_y = sy2
endif

v(1)=(roll_x-cr_x)/dtime
v(2)=(roll_y-cr_y)/dtime

endif
call apply_fixed_bc(nsurf,2,v(1:3))
endif

c if called for centre roll right

if(nsurf.eq.ivarvalue(idx_cr_right)) then
c curernt position of roll
call sfkinbd(ivarvalue(idx_cr_right),xmin,xmax,0)
c Calculate roll center coordinates
cr_x=half*(xmax(1) + xmin(1)) !x
cr_y=half*(xmax(2) + xmin(2)) !y
cr_r=half*(xmax(1) - xmin(1)) !radius

arm_radius=rvarvalue(idx_rarm_rad)

c Get the position of two circle intersection

call twocirc_intersec(0.0d0,wp_y,(wp_r+cr_r),
& cr_r_armx,cr_r_army,arm_radius,
& sx1,sy1,sx2,sy2,nintersec)
if (nintersec.eq.2) then
c rightmost intersection point
if (sx1.lt.sx2) then
roll_x = sx1
roll_y = sy1
else

145
2023.2 User defined kinematic Example

roll_x = sx2
roll_y = sy2
endif

v(1)=(roll_x-cr_x)/dtime
v(2)=(roll_y-cr_y)/dtime
endif
call apply_fixed_bc(nsurf,3,v(1:3))
endif
endif

c=======================================================================
elseif(iflag.eq.2) then
c-----------------------------------------------------------------------
c iflag = 2 End of increment call
c-----------------------------------------------------------------------
c calculate ring diameter and update for THS plot
call sfkinbd(1,xmin,xmax,0)
ring_height=xmax(3)-xmin(3)
if(calc_ring_ypos(z_measure,yy,1)) then !ring profile only along
x=0,z=mid ring
ring_dia=yy(4)-yy(1)

if(inc.eq.0) then
xkinbc(1)=ring_dia ! Outer ring diameter along y axis
xkinbc(2)=ring_dia
else
xkinbc(2)=xkinbc(1) !old diameter
xkinbc(1)=ring_dia
endif

dia_inner=yy(3)-yy(2)
r_thick=half*(yy(4)+yy(2)-yy(3)-yy(1)) !average of left and right
if(inc.eq.0) then
write(kou,8002) "initial ring dimension"
write(kou,8003) ring_dia,dia_inner,r_thick,ring_height
write(kou,8002) "desired ring dimension"
write(kou,8003) rvarvalue(idx_od),rvarvalue(idx_id),
& rvarvalue(idx_wt),rvarvalue(idx_ht)
write(kou,8004)
else
write(kou,8002) "current ring dimension"
write(kou,8003) ring_dia,dia_inner,r_thick,ring_height

if(ring_dia.ge.rvarvalue(idx_od)) then !only outer diameter

critaria is considered
write(kou,8005)
write(iprtl,8005)
iotnum=inc+1 !simulation will stop
endif
endif
endif
call sfukin2_restart(iflag)
call sync_ddm_domain()
c check for possible error and quit
if(kinerror.gt.0) call quit(13)
endif

146
2023.2 User defined kinematic Example

return
7001 format(/13x,'*** error enter non-zero axial roll angle',es21.13)
8001 format(13x,a,' (x,y) :',2es21.13)

8002 format(/,/10x,41('-'),/12x,a,x,'(m)',/10x,41('-'))
8003 format(12x,'outer diameter :',es21.13,/,
& 12x,'inner diameter :',es21.13,/,
& 12x,'ring thickness :',es21.13,/,
& 12x,' ring height :',es21.13)

8004 format(/,/13x,'simulation will end if the above',

& ' outer diameter is reached' )

8005 format(/,/13x,'desired outer dimeter reached.',

& ' simulation will end')

8007 format(/13x,'*** warning the provided measuring point for',

& ' ring diameter at z=',es12.5,' is outside the ring ','(',
& es12.5,',',es12.5,')',
& /13x, 'average ring dimension in z direction will be used' )

8008 format(13x,'ring diameter will me measured at z =',es12.5,

& ' at x = 0 plane along y axis',/ )
end

subroutine sfukin2_restart(iflag)
use sfukin2_module
c-------------------------------------------------------
c function to create kin.res file for restart/prestate
c-------------------------------------------------------
c iflag = 0 to read back if the job is restart
c = 2 to write information on disk for prestate/restart
c
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
#include "op1.cmn"
#include "jname.cmn"
#include "sfkin.cmn"
#include "cdominfo.cmn"
#include "concom.cmn"
#include "marc_prestt.cmn"
#include "iautcr.cmn"

integer iflag,iunit
integer n_int,n_real,n_char
character*(MAX_CARD_DIM_2) filen
logical ddmjob,begin,restart

c unit number
iunit=2021

ddmjob=.false.
begin=.false.
restart=.false.

147
2023.2 User defined kinematic Example

if(inc.eq.0) begin=.true.
if(nprocd.gt.0) ddmjob=.true.
if(iprestt_data(23).gt.0) restart=.true.

c create file name, donot change the name for GUI support
if(ddmjob)then
filen=trim(jidnam_base)//'_kin.res'
else
filen=trim(jidnam)//'_kin.res'
endif

c If restart job
if(iflag.eq.0) then
if(restart) then
c open and read _kin.res
call openfile(iunit,filen,1,5)
call sfkin2res (iunit,0,0,0,0,1,1) !read

c ------------COPY BACK----
c copy back integers
iotnum=i2res(1)

c copy back real/doubles

cr_l_armx = x2res(1)
cr_l_army = x2res(2)
cr_r_armx = x2res(3)
cr_r_army = x2res(4)
tang_vel_mr = x2res(5)
term1 = x2res(6)
yoffset_ini = x2res(7)
xkinbc(1) = x2res(8)
xkinbc(2) = x2res(9)
c -------------------------

c close _kin.res file

call closefile(iunit,filen,0,2)
endif

c if called at end of increment

else if(iflag.eq.2) then

if(iparent.eq.0) then
c ------------SAVE important variables----
c saving integers
i2res(1)=iotnum

c saving real/doubles
x2res(1)=cr_l_armx
x2res(2)=cr_l_army
x2res(3)=cr_r_armx
x2res(4)=cr_r_army
x2res(5)=tang_vel_mr
x2res(6)=term1
x2res(7)=yoffset_ini
x2res(8)=xkinbc(1)
x2res(9)=xkinbc(2)
c update how many integers and real are saved

148
2023.2 User defined kinematic Example

n_int=1
n_real=9
n_char=0
c -------------------------
if(begin .and. .not.restart) then
call openfile(iunit,filen,1,5)
call sfkin2res (iunit,iunit+1,n_int,n_real,n_char,1,2) !write !
change the numbers
else
call openfile(iunit,filen,3,5)
call sfkin2res (iunit,iunit+1,n_int,n_real,n_char,1,3) !append
endif
call flush(iunit)
call closefile(iunit,filen,0,2)
endif !iparent
endif !if(iflag.eq.0) then
return
end

149
Scientific Tutorial
2023.2

3 User defined material

2023.2 User defined material User-defined material law

Educational Objectives

Uses of user subroutine in Simufact Forming for dynamic material properties

Prerequisites

FORTRAN programing knowledge, general programing knowledge

3.1. User-defined material law

Through user subroutines, it is possible to define custom material laws through analytical equations in Simufact
Forming. The flow stress of the material as a function of plastic strain, plastic strain rate and temperature can be de-
fined and calculated during the simulation. This can be achieved through the usage of the subroutines au_yiel() and
user_fstress(). In the case of au_yiel(), all the parameters have to be defined within the subroutine. This routine can be
used if a simple flow stress equation is used for all temperatures. The usage of the subroutine user_fstress() provides
more flexibility and greater process control for complex non continuous temperature dependence of the flow stress. It
is common to have different values of flow stress at different temperatures (cold and hot forming), therefore the ma-
terial may have different constitutive equations for different temperatures. The flow stress of a material is a complex
function of parameters such as strain, strain rate and temperature. In such cases using the subroutine user_fstress()
helps to specify different flow stress equations with different parameters. This manual will focus on the use of the
subroutine user_fstress().

To include the material law in the simulation, two steps have to be performed. First, various parameters for the material
equation have to be determined in the material database through Simufact Material. Second, the material equations
based on these parameters have to be described in the user subroutine user_fstress(). When the parameters are defined
in the material information through Simufact Material, it is important to provide a model number and validity range in
terms of strain, strain rate and temperature for which the model is valid. For example, if the parameter set is assigned
with the model number 1 in Simufact Material, then the solver assigns the parameter lists with the number 101 (model
numbers 1-100 are reserved internally) and will pass these parameters along with the model number 101 (in this case)
through the call arguments in the user_fstress() function. In addition to the parameter list, current strain and strain rate
are also included in the argument list. In case that the current strain or strain rate lies beyond the defined range for the
model, the solver will provide the subroutine with only the upper or lower range value that is nearest to the actual value.

It must be noted that the solver uses temperature ranges to determine the parameter sets that should be passed through
the function arguments. If two different flow stress models are used for two different temperature ranges and there
is discontinuity between the temperature ranges, then the flow stress in the region of discontinuity will be linearly
interpolated from the two end point values of the two models. For example, if equation 1 is valid till temperature
T1_max and equation 2 is valid from temperature T2_min, then the flow stress at a temperature anywhere in between
T1_max and T2_min is obtained by linear interpolation of the flow stress values from equation 1 at T1_max and
equation 2 at T2_min, as shown in the figure below. If the material equations have an overlapping temperature range,
the system will give an error.

151
2023.2 User defined material Material information

Figure 3.1. Flow curves for different temperature ranges and linear interpolation

3.1.1. Material information

A material database for the flow stress model can easily be created through Simufact Material which then can easily
be imported in Simufact Forming. The user has the option to create a custom material from scratch by unlocking the
database. However it is easier to make a copy of the nearest matching material from the existing database and then
change the relevant data. The later process is recommended because less effort is required to add or change the material
properties, other information remain unchanged. It is not recommended to edit and overwrite the already existing
material properties. So, first make a copy and then make the changes. A material properties window looks like below:

Figure 3.2. Flow curve information created through table input

152
2023.2 User defined material Material information

3.1.1.1. New material database

In this section, the steps to create a new material database from an already existing one are described and an example
is given.

In Simufact Forming, first import a material in the inventory window. The imported material should be the nearest
matching to that being used for the simulation process so that only selected few information has to be altered. In this
example, we import stainless steel with DIN number 1.4571 which can be used for cold as well as for hot forming.
The figure below shows the inventory window with a material that can be duplicated through right mouse click and
can be renamed afterwards.

Figure 3.3. Inventory window showing process of duplicating the material file

The properties window of the new material can be opened by a double mouse click on it. The window shows various
material properties that can be modified. For the flow curve menu on the left, it shows the existing data which can be
removed by clicking on the minus "-" button on the right.

153
2023.2 User defined material Material information

Figure 3.4. Removing the existing flow curve data

After removing the existing flow curve data, choose the Flow curve approach Analytical in the drop down menu,
as marked with 1 in the figure below. Then choose User-defined in the Plasticity model option and confirm the
choice by pressing the green "+" button, as shown by number 2 and 3 respectively in the figure below.

154
2023.2 User defined material Material information

Figure 3.5. Adding a user-defined material plasticity model

In the following dialog window various parameters can be entered to describe the flow curves. Press the edit button
to open the dialog window:

155
2023.2 User defined material Material information

Figure 3.6. Open the dialog window to edit the flow curve data

The dialog window for the parameter definition looks like in the figure below. The values for various parameters can
either have a constant value or can be entered through a table.

Figure 3.7. Dialog window to enter equation parameters

156
2023.2 User defined material User subroutine

To enter the values in the above dialog window follow the following steps:

1. Press the green '+' button beside the temperature range option to activate the input fields.

2. Select the appropriate unit of temperature and stress from the drop down menu

3. First enter the temperature range for which the parameters specified in this window are to be used in the calculations.
Next enter the minimum and maximum values of strain and strain rate for which the model will be valid. In addition
total of maximum 16 different parameter values can be entered but not all have to be used. The values of these
parameters can be a constant or a function of time, temperature and strain rate. The relation has to be entered
through a table input. The solver will give an error if two equations have an overlapping temperature range, so
please avoid it.

4. Enter the model number. The solver will use the subroutine argument itype1 to pass the variable array set depending
upon the model number. itype1 is 100 plus the model number entered here. Because of this it is important to assign
a distinct model number.

5. Save the entries by clicking the apply button. To add parameters for another model, repeat the above steps from 1
- 4. Make sure that different model numbers are used for different parameter sets.

Different temperature ranges can share either the same model number (i.e. the same equation defined in user_fstress.f)
using different parameters or different model number. Both is supported and both may be reasonable depending on
the behavior of your material and your design of the user subroutine. Take care to define the material accordingly.

The new created material database can be saved in the Simufact Material database for future use through the "To
Library" option in the inventory window.

Figure 3.8. Save the created material file for future use

3.1.2. User subroutine

The function user_fstress() will be used to define the analytical material law. The subroutine is called for each inte-
gration point of the element where the flow stress of the material has to be calculated. The syntax of the subroutine
is as follows:

subroutine user_fstress(ilst,ityp1,t1,ep1,er1,p1,fstress1)
c
c subroutine to compute flow stress based on user equation
c ilst(1) - material id

157
2023.2 User defined material Example

c ilst(2) - user element number

c ilst(3) - integration point
c ilst(4) - layer id
c ilst(5) - internal layer number
c ityp1 - equation type (>100)
c Note
c This is the "Model number" defined in the GUI + 100
c E.g. If "2" is defined in the GUI, ityp1 is 102
c t1 - temperature
c ep1 - plastic strain
c er1 - plastic strain rate
c p1 - material parameters from input
c fstress1 - flow stress
c
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif
real*8 t1,ep1,er1,p1,fstress1
integer ilst,ityp1
dimension p1(*),ilst(*)
cccccccccccccccccccccccccc

c USER DEFINED EQUATIONS

cccccccccccccccccccccccccc
return
end

Through the above subroutine, the various parameters defined in the material information are passed to the function. In
the argument, the solver passes the model number through the variable ityp1 and parameter values through arrays p1(),
like p1(1), p1(2) and so on for parameter 1 and parameter 2 respectively. At the same time the current temperature,
plastic strain and strain rate are available through the variables t1, ep1 and er1 respectively. The solver will use the
temperature field to determine which equation number and corresponding sets of parameters to be made available
through the function arguments.

3.1.3. Example
The subroutine shown below calculates the flow stress through a simple power law. Two equations are defined to
calculate the flow stress. The equation number is chosen on the basis of the temperature. At low temperature the first
equation will be used whereas the second one is used for high temperatures.

subroutine user_fstress(ilst,ityp1,t1,ep1,er1,p1,fstress1)
c
c subroutine to compute flow stress based on user equation
c ilst(1) - material id
c ilst(2) - user element number
c ilst(3) - integration point
c ilst(4) - layer id
c ilst(5) - internal layer number
c ityp1 - equation type (>100)
c Note
c This is the "Model number" defined in the GUI + 100
c E.g. If "2" is defined in the GUI, ityp1 is 102
c t1 - temperature
c ep1 - plastic strain
c er1 - plastic strain rate

158
2023.2 User defined material Example

c p1 - material parameters from input

c fstress1 - flow stress
c
#ifdef _IMPLICITNONE
implicit none
#else
implicit logical (a-z)
#endif

real*8 t1,ep1,er1,p1,fstress1
integer ilst,ityp1
dimension p1(*),ilst(*)

c------Equation 1
if(ityp1.eq.101) then
fstress1=p1(1)*er1**p1(4)
if(fstress1.le.p1(2))fstress1=p1(2)
endif

c-------Equation 2
if(ityp1.eq.102) then
fstress1=p1(1)*er1**p1(4) * ep1**p1(3)
if(fstress1.le.p1(2))fstress1=p1(2)
endif
return
end

In this example, the same simulation is set up with different temperature objects (500°C respectively 1000°C) for the
workpiece and the dies. Thus different equations are used to get the flow curve for the two cases. Two parameter sets
are defined in the material file for two temperature ranges as shown below:

Figure 3.9. Example - Parameters for temperature range 1 (300°C-600°C)

159
2023.2 User defined material Example

Figure 3.10. Example - Parameters for temperature range 2 (600°C-1200°C)

Comparing the z-forces of the upper dies for the two cases highlights the use of different flow stresses due to different
simulation temperatures.

Figure 3.11. Example - z-forces of the upper dies

160
Hexagon is a global leader in digital reality solutions, combining sensor,
software and autonomous technologies. We are putting data to work
to boost efficiency, productivity, quality and safety across industrial,
manufacturing, infrastructure, public sector, and mobility applications.

Our technologies are shaping production and people-related ecosystems

to become increasingly connected and autonomous – ensuring a scalable,
sustainable future.

Hexagon’s Manufacturing Intelligence division provides solutions

that use data from design and engineering, production and metrology
to make manufacturing smarter.

Learn more about Hexagon (Nasdaq Stockholm: HEXA B) at hexagon.com

and follow us @HexagonAB.

Copyright © 2023 Hexagon AB and/or its subsidiaries. All rights reserved. Hexagon, the Hexagon logo, and
other logos, product and service names of Hexagon and its subsidiaries are trademarks or registered
trademarks of Hexagon AB and/or its subsidiaries in the United States and/or other countries. All other
trademarks belong to their respective owners.