FEAP User Forum
FEAP => Programming => Topic started by: R2197 on August 08, 2023, 02:04:22 AM
-
Hi there I have a quick question.
I want to add a subroutine in my user defined element. However, this should be outside the time loop as I want to perform some computation just once and then call the output in the time loop iterations. Now, as I understand the time steps and loops are governed by batch files. So, I am not exactly sure where to put my subroutine so that it runs only once.
Thanks
-
You can try putting it as part of isw=14. That is called once after mesh input.
-
Dear Prof. Taylor
Thanks for your reply. But, I am not quite sure how to execute this. Whether I should add this in the material or element files and if something like an if condition below enough or not?
if (isw.eq.14) then
call subroutine(a,b,c)
end if
Are there any example codes that I can have a look at?
Also maybe it is better to explain the particular problem. I am trying to implement a FFT based solver using the FFTW3 library. The implementation works but now I want to separate the 'plan' creation and take it outside the time loops. So, I intend to create the plan only once and call it for all subsequent iterations. This would save me some computational effort.
-
isw.eq.14 is passed to user elements at the start of a computation. it occurs sometime after the mesh data has been read. typically it is used to initialize history variables and other such items that only have to be done once during analysis.
-
Describe what data the plan creates. Do you make the plan in the element routine now?
-
Dear Prof. Govindjee and Prof. Taylor
The plan subroutine creates the best possible way to ocmpute the FFT of the given array. A few things to consider while creating the plans are for eg. size of the array, whether to compute the FFT of N sized array or break it down into P smaller problems of N/P size each, the way to combine the results if the solution is determined by breaking down the large N array, etc. All of this information is in data structures stored as pointers whose location is the value of the 'plan' variable.
I was using a subroutine in the material file. And this subroutine was calling another routine for calculating FFTs and the planner was also a part of this. Following is the flow:
Material File
-call FFT()
-call Planner()
But the issue with this type of implementation is that it calls the planner every time which consumes a lot of time and memory. I want to do it only once before the time loop.
So, as suggetsed before that isw.14 might work, I am now looking into how to create a umacrX file with a prompt, let's say 'PLANNER' which can be called with the other BATCH commands.
-
Yes, using a umacrXX is another way to do this. Then you can call the routine whenever you want from the command prompt.
-
Dear Prof. Govindjee
Thanks for your comment.
But I was still wondering how will this plan variable then be called by the FFT() subroutine. Can it just call it in directly or do I have to define initialise or do something there as well?
-
From your description you want to call the plan as part off the material model routine.
My suggestion is to use isw=14 to do this. You merely use your element routine up to where the material model is done and do the plan. You do not need to do the rest of the computation but could as the element stiffness and residual will not be used in isw=14. You will need isw passed to your material routine so you can either compute or use the plan.
If you use a umacrxx you would need to call the element routine yourself. Of course this could be done using the routine 'form' but then you have to code the element routine to use whatever isw you choose. So no advantage over the isw=14.
-
Dear Prof. Taylor
I was not exactly sure which is better and so thought about the UMACRXX. But I can also try isw.eq.14. The thing is I am not sure how and where to implement it in the element file. Should I do it around line 157 in this file?
The function matlib_AS just calls the material file according to the ID and its relevant data.
-
You have a large IF-ELSEIF-END structure:
if(isw.eq.1) then
...
elseif(isw.eq.3.or.isw.eq.4.or.isw.eq.8) then
...
end if ! isw: 1,3,4,8
Place your elseif(isw.eq.14) then
...
code just before the end if ! isw: 1,3,4,8
-
It would help if you could describe what needs to be known for you to create a plan.
You placed the commented out plan inside an isw=1 part, so that could never work with the "c" removed.
You indicated originally that you had a working module with a plan, but it created the plan every iteration. Can you show that program module? From that we should be able to help more.
-
Dear Prof. Taylor
The commented out plan in isw.eq.1 was my attempt before asking in the forum. I am attaching another file which has the relevant code for the Fourier Transforms including plan creation(lines 84-91) and execution(lines 148,159 etc.). The function endings R2C and C2R are from real to complex transforms and vice versa. What I am trying to do now is take these plan creation routines out of this subroutine and perform them only once as I do not change the size of the arrays which is given from the material file at the moment.
To create the plan, one has to provide the size of the transform, the input and output arrays, and flags indicating the desired behavior of the plan (e.g., FFTW_FORWARD or FFTW_BACKWARD for the transform direction). Then the planner function searches the set of plans and provides an optimum plan considering also the plan creation flags(FFTW_EXECUTE, FFTW_MEASURE, etc.). And quoting from "Implementing FFTs in Practice" by S. G. Johnson and Matteo Frigo, the theoretical basis behind these plans is, "Roughly speaking, to solve a general DFT problem, one must perform three tasks. First, one must reduce a problem of arbitrary vector rank to a set of loops nested around a problem of vector rank 0, i.e., a single (possibly multi-dimensional) DFT. Second, one must reduce the multi-dimensional DFT to a sequence of of rank-1 problems, i.e., one-dimensional DFTs; for simplicity, however, we do not consider multi-dimensional DFTs below. Third, one must solve the rank-1, vector rank-0 problem by means of some DFT algorithm such as Cooley-Tukey. These three steps need not be executed in the stated order, however, and in fact, almost every permutation and interleaving of these three steps leads to a correct DFT plan."
-
Is it possible to add a "switch" such that after the first call to the plan the switch is changed.
Thus in the first call the plan is executed and saved. In succeeding calls the plan is read from the saved file. The execution is always performed.
One such switch is a logical variable. Initially set to "true". Thus, the code is
If (switch) then
do plan and save
switch = .false.
else
read plan
endif
You can initialize the switch during input of material set data possibly. Put in module or common block to pass to the plan code you posted.