I was doing stress analysis for an aeroplane wing and I wanted to apply a parabolic distributed load along the length of the wing… Can someone explain how to implement such a non uniform distribution?
See Section 8.4 of the ccx manual:
8.4 User-defined loading: These routines are made available to define nonuniform distributed loading. The user can define the loading in each integration point separately as a function of position, time etc.
8.4.1 Concentrated flux (cflux.f)
This subroutine is used for user-defined concentrated heat flux, characterized by the parameter USER on the *CFLUX card. For the header and variable description, the reader is referred to cflux.f in the source code.
8.4.2 Concentrated load (cload.f)
This subroutine is used for user-defined concentrated load, characterized by the parameter USER on the *CLOAD card.
When using a *DLOAD card, you can actually put a pressure on individual element faces.
So that’s what I used when I had to implement a non-even distributed load. (in my case, an IEC 60601 human body mass distribution for a medical device.)
To make this work, I had cgx write out the load surface as a *SURFACE card next to the nodes and elements (all.msh).
A Python script was used to determine all the facets that carried part of the load based on their position and wrote them out as a *DLOAD for the individual facets.
The facets that carry part of the load are indicated in blue below.
(Note that this doesn’t show the magnitude of the individual pressures, which can be different for each facet.)
In your case you could calculate the center of each facet (by averaging the coordinates of the corner nodes) and assign a pressure to that facet based on the location of the center.
Thanks for your help, I really like that way but I still don’t get how did you assign different loads to different elements? How did you tell python the load that should be carried by each element
Tbh, I don’t use Fortran so I dont know how to type such subroutines, do you have an example so I can focus on it and try to figure out how is it written?
It is simply specified as a DLOAD, for example;
*STEP
*STATIC
*DLOAD
7349,P2,3064.83
7350,P2,5994.72
...
8120,P2,51573.05
*NODE FILE
U
*EL FILE
ZZS,ME
*END STEP
For simplicity I stick to SI units in my FEA.
So, a pressure of 3064.83 Pa acts on face 2 of element 7349, et cetera.
How to calculate the pressure to use depends on the situation.
Say for example one has an even pressure on a circular surface. For the element faces that are completely inside the circle, the pressure to use is simply the pressure on the circular surface. For elements that are partially inside the circle, you multiply the pressure by the fraction of the face area that is inside the circle. That fraction is determined in a simplified manner by simply mapping the amount of nodes from the facet to a surface fraction, using a Python dictionary. This is for a C3D20 element, that has 8 nodes in every facet.
# The key is the number of nodes that fall inside of a geometry.
# The value is the average area of a facet that is then occupied.
SURFACE_FRACTION_C3D20 = {
1: 1 / 16,
2: 1 / 8,
3: 1 / 4,
4: 1 / 2,
5: 3 / 4,
6: 13 / 16,
7: 15 / 16,
8: 1,
}
when subroutine is need to be eliminated, probably an approximation by method of averaging and finer mesh still can be acceptable in results.
parabolic pressure can be approximate by large piecewise of linear values. Maximum and minimum of pressure values at each edge being averaged instead of actual values, structured mesh is required to manage easily in solid face of group definition.
another problem may raise when pressure at longitudinal direction is not consistent, the task become cumbersome. Imported from OpenFoam result directly can simplify if available.
Hi Aly,
I don’t have programming skills so my suggestion is more pedestrian, but it works.
Let’s consider the wing centered and aligned along the X axis.
Apply a unit traction load on the wind (with Prepomax or whatever software you use,…).
Let the software generate the inp for you. That will help you to know the weight of each node when the load is distributed and identify which nodes you need to consider.
Edit the inp file with Excel or your favorite Spreadsheet. (Or python,…)
Now you need to multiply each node under the CLOAD card shown on the inp by the new parabolic pressure profile ( Ax^2+B*x+C) where x is the x coordinate of each node under consideration.
For example, starting from this inp file generated by Mecway in this case. I have considered a traction of 1KPa over the surface of interest (C3D20R):
*CLOAD
14,3,-0.01666666666667
16,3,-0.01666666666667
15,3,-0.01666666666667
13,3,-0.01666666666667
97,3,0.03333333333333
88,3,0.06666666666667
98,3,0.03333333333333
.
You can clearly see how the weight of each node is different when the load is distributed. Then compute the new CLOAD in Excel like this:
*CLOAD
14,3,-0.01666666666667*( A*X_14^2+B* X_14+C)
16,3,-0.01666666666667*( A*X_16^2+B* X_16+C)
15,3,-0.01666666666667*( A*X_15^2+B* X_15+C)
13,3,-0.01666666666667*( A*X_13^2+B* X_13+C)
97,3,0.03333333333333*( A*X_97^2+B* X_97+C)
88,3,0.06666666666667*( A*X_88^2+B* X_88+C)
98,3,0.03333333333333*( A*X_98^2+B* X_98+C)
.
.
.
A, B depends on your desired parabolic shape being x^2 the simplest one.
C parameter is used to scale the whole load. Check overall reaction to adjust the value.
This is the result:
2 Likes
If no experience with Fortran or subroutines, then don’t waste your time trying that method. I believe your best bet is the solution by @Disla below.
That’s a great simple conceptual solution.
One important detail; the lift on an airplane wing is largest in the centre and tapers off to the ends. So it’s more like an inverted parabola.
I really appreciate your replies, and yes Disla’s solution is really nice but one last question how can you generate that (node,force value) table using cgx? what I know is we define a set of nodes (for example “load”) then we type in a text file *CLOAD followed by load,3,1… so how can we show up the node number with the force so I can multiply the unit load by a function?
typically a 1D model is used to apply external loads and generate internal forces (shear force, bending and torsional moment mainly). This model is made up of beams with same EI and GJ as the 3D wing (as accurate as possible, will evolve during design cycle) and place in the elastic centre of the wing. As is a 1D model 1 single coordinate define the load distribution son you only need to know nodal coordinates along the wing span and apply the magnitude according to your distribution.
Remember that aircraft wing loading is elliptic, not parabolic, with some reduction due to fuselage interference.
In the 3D model internal forces generated are later applied using compliant links (usually attached at rib locations).
sorry, previously i miss to interpreted about parabolic pressure distribution at longitudinal direction. I guessed before, something like wind load on bridge a non uniformity of pressure apply in section.
(source: SOFiSTiK AG & Cyclone Fluid Dynamics BV, 2012)
1 Like
I would like to thank you all for explaining and giving different solutions. Among other softwares i used, this forum have always been the best
3 Likes
Cantilever beam under parabolic load.
| Cantilever under parabolic load | |||||||
|---|---|---|---|---|---|---|---|
| L | 7.39 | [m] | |||||
| h | 0.20 | [m] | |||||
| P0 | 24,825 | [N] | Result | 122,305 [N] | - 339,999 [Nm] | - 0.280 [m] | |
| E | 210,000,000,000 | [Pa] | Theory | - 122,305 [N] | 338,936 [Nm] | 0.279 [m] | |
| I | 0.000066667 | [m4] | Deviation | 0.00% | 0.31% | 0.38% |
And this is the reference in case someone wants to try the procedure. (Edit inp + some Excel)
