Map 3D Soil Variabilities into STKO soil model

sebi_magi
Posts: 21
Joined: Mon May 30, 2022 5:50 pm

Re: Map 3D Soil Variabilities into STKO soil model

Post by sebi_magi » Wed Jul 06, 2022 6:43 am

Hello there,

I attach my python code that does what I want. It works well but it takes around 30 minutes when I evaluate over the whole geometry.

Code: Select all

from PyMpc import *
from PyMpc import MpcOdbVirtualResult as vr
import numpy as np
import os
#import pandas
# clear terminal
App.clearTerminal()

# get document
doc = App.postDocument()
# get first database
if len(doc.databases) == 0:
		raise Exception("You need a database with ID = 1 for this script")
def get_Nodal_disp(finput,fdispx,fdispy,fdispz,fnodes,nodresult='Acceleration'):
	db=doc.getDatabase(2)
	
	nodes_target=np.loadtxt(finput,dtype=int)
	print('n',nodes_target)
	displacement = db.getNodalResult(nodresult, match = MpcOdb.Contains)
	# create evaluation options
	# here we want to extract data for all steps of the last stage
	all_stages = db.getStageIDs() # get all model stages
	last_stage = all_stages[-1] # get the last stage
	all_steps = db.getStepIDs(last_stage) # get all steps of the last stage
	all_times = db.getStepTimes(last_stage) # get all step times of the last stage
	# initialize the evaluation options with the last stage
	opt = MpcOdbVirtualResultEvaluationOptions()
	opt.stage = last_stage
	num_steps = len(all_steps)
	## create node array with disp for all nodes, first row is time
	arr_disp_x=np.zeros((num_steps,len(nodes_target)+1))
	arr_disp_y=np.zeros((num_steps,len(nodes_target)+1))
	arr_disp_z=np.zeros((num_steps,len(nodes_target)+1))

	## create node array with node id and coordinates
	arr_nodes=np.zeros((len(nodes_target),4))
	# flag for updating arr_nodes only at first time step
	create_node_file=True
	for step_counter in range(num_steps):
		# get step id and time
		step_id = all_steps[step_counter]
		step_time = all_times[step_counter]
	
		# write something...
		# and process all application events to avoid GUI freezing...
		print('Evaluate results at step {} t = {}'.format(step_id, step_time))
		App.processEvents()
	
		# put the current step id in the evaluation options
		opt.step = step_id
	
		# evaluate all the results at the current step
		displacement_field = displacement.evaluate(opt)
	
		# This is how you extract data from a nodal result at one or multiple nodes
		# and compute the relative drift accessing the mesh data
		# this is a nodal field. The row is the nodal id
		all_disps = []
		for n,nod in enumerate(nodes_target):
			## get nodal id
			i_node_id=nodes_target[n]
			# not sure what this does but works
			i_row=MpcOdbResultField.node(int(nod-1))
			Ux=displacement_field[i_row,0]
			Uy=displacement_field[i_row,1]
			Uz=displacement_field[i_row,2]

			## update time step
			arr_disp_x[step_counter,0]=step_time
			arr_disp_y[step_counter,0]=step_time
			arr_disp_z[step_counter,0]=step_time

			## update disp
			arr_disp_x[step_counter,n+1]=Ux
			arr_disp_y[step_counter,n+1]=Uy
			arr_disp_z[step_counter,n+1]=Uz

			## if the flag for arr node coords is True, then create nodedata
			if create_node_file:
				node_data=displacement_field.mesh.nodes.values()[nod-1]
				#print('bu',nod,n)
				#print('x',node_data)
				nodepos=np.array([nod,node_data.x,node_data.y,node_data.z])
				arr_nodes[n,:]=nodepos
		## turn off the flag!
		create_node_file=False
	# create a string with the node disp header
	string=[str(x) for x in nodes_target]
	string.insert(0,'time')
	header_data=';'.join(string)
	## create a string with the node coords header
	header_nodes='id;x;y;z;'
	#print(arr_disp_x.shape)
	#print(arr_nodes)
	## format for node_coords output file
	format_nodes=['%i','%10.5f','%10.5f','%10.5f']

	np.savetxt(fdispx,arr_disp_x,fmt='%10.5f',header=header_data,comments='',delimiter=';')
	np.savetxt(fdispy,arr_disp_y,fmt='%10.5f',header=header_data,comments='',delimiter=';')
	np.savetxt(fdispz,arr_disp_z,fmt='%10.5f',header=header_data,comments='',delimiter=';')

	np.savetxt(fnodes,arr_nodes,fmt=format_nodes,header=header_nodes,comments='',delimiter=';')

lay='surface'
sta='AP01'
disp_type='V'
finput='MAULE_NODES_'+lay.upper()+'.txt'
fdispx=sta+'_'+disp_type+'x_'+lay+'.csv'
fdispy=sta+'_'+disp_type+'y_'+lay+'.csv'
fdispz=sta+'_'+disp_type+'z_'+lay+'.csv'
fnodes=sta+'_nodes_'+lay+'.csv'
get_Nodal_disp(finput,fdispx,fdispy,fdispz,fnodes,nodresult='Velocity')
print('done')


I'm also having issues also with modelling cross sections from the 3D Model. I notice that this slice of the 3D model won't be alligned with the global coordinate axis. When I run the model to map the soil properties, I get an error like this. Perhaps there is an issue with the normal vectors of the absorbing boundaries
cross_section_with_error.png
cross_section_with_error.png (232.33 KiB) Viewed 3466 times
Kind regards.

sebi_magi
Posts: 21
Joined: Mon May 30, 2022 5:50 pm

Re: Map 3D Soil Variabilities into STKO soil model

Post by sebi_magi » Mon Jul 11, 2022 5:01 am

Hello,

I managed to solve most of the problems I had. It was much easier to rotate the cross-sections outside STKO and then import the geometry. This is a 2D Shear Wave velocity profile plot I made. Your random material wrapper and plotting tools are wonderful. I really appreciate your work.
perfil_2d_1.png
perfil_2d_1.png (25.88 KiB) Viewed 3446 times
Then, I input a very large earthquake to my model. The input (red), bottom (dark) and surface (blue) accelerations are shown in the next figure. I averaged bottom and surface nodal accelerations to get the red an dark Time Series.


The fact that the input is no the same than the nodal bottom acceleration bugs me. I think it is because multiple waves are propagating and reflecting through the soil, so we cannot expect the ideal case of the ricker wavelet experiment propagating through homogeneous media (I run a ricker wavelet experiment and get coherent results). Do you agree?
amplificacion.png
amplificacion.png (37.85 KiB) Viewed 3446 times
After this I want to use an elastic-plastic constitutive model for the 2D profiles. Do you think it will be problematic with the ASDA Absorbing boundary?

Lastly, does the element thickness plays a role in the 2D models? I dont know what varying the element thickness does.

Kind regards.

STKO Team
Posts: 3066
Joined: Tue Oct 29, 2019 8:45 am

Re: Map 3D Soil Variabilities into STKO soil model

Post by STKO Team » Mon Jul 11, 2022 8:08 am

You did a great job solving those issues by yourself!
managed to solve most of the problems I had. It was much easier to rotate the cross-sections outside STKO and then import the geometry. This is a 2D Shear Wave velocity profile plot I made. Your random material wrapper and plotting tools are wonderful. I really appreciate your work.
This is true, because in OpenSees a 2D analysis can ONLY BE in the XY plane!
The fact that the input is no the same than the nodal bottom acceleration bugs me. I think it is because multiple waves are propagating and reflecting through the soil, so we cannot expect the ideal case of the ricker wavelet experiment propagating through homogeneous media (I run a ricker wavelet experiment and get coherent results). Do you agree?
Exactly, you cannot see the exact same input, because there are also the waves reflected from the top free surfaces that sum up! With the wavelet it's easy because by the time the wave gets reflected, the input at the base is already done
After this I want to use an elastic-plastic constitutive model for the 2D profiles. Do you think it will be problematic with the ASDA Absorbing boundary?
Do it! You will find an issue when the model goes in plastic state, because the G modulus of the absorbing boundary remains constant. It is possible now to update it step by step with a tcl script. Try and we will guide you in how to do it.
FYI, we are implementing a new feature in the absorbing boundary: You will be able to input a material instead of G,v,rho constants! It will be available soon

Lastly, does the element thickness plays a role in the 2D models? I dont know what varying the element thickness does.
No, it will just affect the magnitude of the reactions!

sebi_magi
Posts: 21
Joined: Mon May 30, 2022 5:50 pm

Re: Map 3D Soil Variabilities into STKO soil model

Post by sebi_magi » Fri Jul 15, 2022 8:28 am

Hello there,

I am having some issues when mapping to the Absorbing boundaries in the Pressure Dependent Multi Yield Model (PDMY).
I want to map the parameter 'RefShearModul' from the material model which is exactly what the absorbing boundaries needs (G). I get the following error:
error_mapping.png
error_mapping.png (103.45 KiB) Viewed 3424 times
Perhaps I am calling it with the wrong name. I attach my .scd and .rmt files.

Kind regards,

Sebastián
Attachments
PDMY_RANDOM.rar
(340.12 KiB) Downloaded 178 times

STKO Team
Posts: 3066
Joined: Tue Oct 29, 2019 8:45 am

Re: Map 3D Soil Variabilities into STKO soil model

Post by STKO Team » Fri Jul 15, 2022 12:54 pm

Unfortunately the PDMY material does not behave like the others... it has a parameter called shearModulus, but it can be only changed, you cannot get it from the parameter. Unfortunately it has been programmed not following the guidelines used for other materials.

However there is a workaround, you can ask for the element material tangent matrix (which coincides with the elastic one at that stage of the analysis) and from it I get G, Bulk, and poisson.

Another issue you may have is that in this way you still cannot get the rho. I also noticed that this material model does not implement the NDMaterial::gerRho() method.. so it will not pass the mass density to the element itself!!! make sure to put the rho also in the quad element
PDMY01_RANDOM_TEST.zip
(254.45 KiB) Downloaded 185 times

sebi_magi
Posts: 21
Joined: Mon May 30, 2022 5:50 pm

Re: Map 3D Soil Variabilities into STKO soil model

Post by sebi_magi » Mon Jul 18, 2022 3:10 am

That was definitely not intuitive :shock: It worked well though, thank you very much for all your answers :)

For the time being I won't try more complex models but rather optimize what I have done, but I have some last questions.

i) When you said I could update G for each time step so absorbing boundaries match soil, a very bad solution I thought is running analysis of 1 time step, then retrieving the updated shear modulus, then continue with the next step. Is there a workaround or this is actually the only way to update G?

ii) See here a Figure from the Paper Liquefaction and Cyclic Softening at Balboa Boulevard
during the 1994 Northridge Earthquake
from Pretell et al., 2021. As you can see, the soil is heteregoneous with patches of sand-like soils and clay-like soils, which are modelled with different constitutive models. For replicating something like this it would be cool if the RandonMaterialWrapper could accept as an input material models! Something like this

1 ElasticIsotropic rho1 G1 v1
2 ElasticIsotropic rho2 G2 v2
3 PDMY rho, refShearModul, refBulkModul, etc..
4 PM4SAND Dr G0 rho, etc..

and then map all the x,y,z with the corresponding material model...

Perhaps this is easier said than done! but it would be a really useful implementation for geotechnical modelling.

Kind regards.
Attachments
figure.png
figure.png (894.05 KiB) Viewed 3398 times

STKO Team
Posts: 3066
Joined: Tue Oct 29, 2019 8:45 am

Re: Map 3D Soil Variabilities into STKO soil model

Post by STKO Team » Mon Jul 18, 2022 8:41 am

i) When you said I could update G for each time step so absorbing boundaries match soil, a very bad solution I thought is running analysis of 1 time step, then retrieving the updated shear modulus, then continue with the next step. Is there a workaround or this is actually the only way to update G?
There is a much easier way.
It's not yet documented because we use it only internally now. But we are planning to make a webinar on this. Basically, if you look the generated Tcl files in analysis_steps.tcl, we do not call ANALYZE $N (where $N is the number of steps). Instead we call ANALYZE 1, in a loop. Before and after each ANALYZE command, we call custom procedures.
If the user wants to do some custom operation either before or after, it can be easily done.
If you are interested, please prepare a working file with the PDMY material with a simple Ricker Wavelet as an input (so that we can do a benchmark). Make it strong enough to introduce a good level of nonlinearity in the PDMY material.
You will see that the more the PDMY goes into nonlinearity, the less accurate the ability of the absorbing boundary to absorb waves (because the G in the PDMY is changing).
Send it here, I will put custom procedures to do the update.
You can start with the file PDMY01_RANDOM_TEST.zip you sent in the previous post, so that we already have the script for coupling the absorbing element to the nearest solid element

it would be cool if the RandonMaterialWrapper could accept as an input material models! Something like this
That's a super good idea! I already opened a new Issue to implement it. It will be very easy to do. In the future, whenever you have requests like this, you can ask it in the forum, but also create a new Issue on our github web: https://github.com/ASDEAsoft/external_s ... TKO/issues

sebi_magi
Posts: 21
Joined: Mon May 30, 2022 5:50 pm

Re: Map 3D Soil Variabilities into STKO soil model

Post by sebi_magi » Sat Jul 23, 2022 4:22 am

Dear STKO Team,

Thank you for taking into account my suggestion :).

Here is the file .scd and .rmt files for a PDMY model with a ricker wavelet. I used a cfactor of 3 so the acceleration input is 3g. I noticed undesired amplifications occur in the boundaries with such a big input, and they are less noticeable by decreasing the accel. input.

Kind regards.
Attachments
RICKER_TEST.rar
(198.69 KiB) Downloaded 222 times

STKO Team
Posts: 3066
Joined: Tue Oct 29, 2019 8:45 am

Re: Map 3D Soil Variabilities into STKO soil model

Post by STKO Team » Mon Aug 01, 2022 1:13 pm

Dear user, we have seen that updating G from Tcl is very difficult and error prone, because the updated G depends on the material model you are using (PDMY01 in your case) and on some internal variables.

We are planning to update the source code of the ASDAbsorbingBoundary in OpenSees such that you can assign a real material model, instead of G,v and rho. In this way you will be able to directly assign the same material you assign to the interior soil domain.

We will let you know as soon as it is implemented.

sebi_magi
Posts: 21
Joined: Mon May 30, 2022 5:50 pm

Re: Map 3D Soil Variabilities into STKO soil model

Post by sebi_magi » Tue Aug 02, 2022 2:02 am

Looking forward to it :) in the meantime I can try other boundary conditions.


Kind regards.

Post Reply