Envelope of accelerations in OpenseesMP with API

Hello! In this section we'll talk about the post processing module and the data visualization.
Post Reply
gabrifio
Posts: 18
Joined: Wed Dec 16, 2020 12:35 pm

Envelope of accelerations in OpenseesMP with API

Post by gabrifio » Tue Apr 18, 2023 1:25 pm

Dear STKO team,

I am trying to use Python API to visualize the envelope of minimum and maximum accelerations for each node of the mesh. I have followed the example of the Getting started Week 12 example, Utils 14 and 15 (load combination and envelope).

In my case I just want to get the envelope of the accelerations (it is a transient (seismic analysis)).

If I understand well, the basic code should be something like this:

Code: Select all

from PyMpc import *
from PyMpc import MpcOdbVirtualResult as vr

App.clearTerminal() 
doc = App.postDocument() 

db = doc.getDatabase(1) 

result_names = ['Acceleration']
result_labels = [
	['Ax','Ay']]


R1 = db.getNodalResult(
	'Acceleration (Volumes',
	match=MpcOdb.StartsWith,
	time=vr.stageEnd(3)
	)
	
print(type(R1))

env_max = vr.max(R1)
env_min = vr.min(R1)

env_max.setDisplayName('(ENVmax) Acceleration')
env_max.setComponentLabels(['Ax','Ay'])

env_min.setDisplayName('(ENVmin) Acceleration')
env_min.setComponentLabels(['Ax','Ay'])

db.addVirtualResult(env_max)
db.addVirtualResult(env_min)
Until line 23 "print(type(R1))" I don't get any errors, then for env_max I get the following error:
Capture.JPG
Capture.JPG (44.99 KiB) Viewed 1489 times
Is this a problem connected with the fact that I am using OpenSeesMP?
How can I modify the code to make it work?

Thank you very much

gabrifio
Posts: 18
Joined: Wed Dec 16, 2020 12:35 pm

Re: Envelope of accelerations in OpenseesMP with API

Post by gabrifio » Tue Apr 18, 2023 1:27 pm

Just for the sake of clarity, I am only interested in the component "Ax".
Thank you again

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

Re: Envelope of accelerations in OpenseesMP with API

Post by STKO Team » Wed Apr 19, 2023 3:03 pm

The error is saying that you called max (and min) with only one argument:

Code: Select all

env_max = vr.max(R1)
Instead in the util you are referring to, we were doing the envelope of 2 results:

Code: Select all

env_max = vr.max(R1, R2)
env_min = vr.min(R1, R2)
But if I understood correctly you need to do the envelope of a single result in time.
Have a look at webinar #20: Intro to STKO v2.0's Python API

gabrifio
Posts: 18
Joined: Wed Dec 16, 2020 12:35 pm

Re: Envelope of accelerations in OpenseesMP with API

Post by gabrifio » Thu Apr 20, 2023 12:15 pm

Dear STKO team,

Thank you I have watched the seminar and used the example test_virtual_result_query_node to change my code.

Now I have modified the code like this. Since I want the results for all the nodes, and I couldn't find a way to get the full number of nodes, I have just defined NN=6502 which is my total number of nodes,

The code seem to run now with no errors, except this one that appears in the terminal: "'MpcOdbVirtualResult' object is not subscriptable"

Could you please help me to make it work? Thank you!

Note that I have created a progress bar to monitor the progress of the calculation since I have a lot of nodes, the dialog bar opens but then it seems that the code is stuck and does not progress. I don't know if this can be useful, in case it can be cancelled.

Code: Select all

import numpy as np
from PyMpc import *
from PyMpc import MpcOdbVirtualResult as vr
from time import sleep
import traceback
from PySide2.QtCore import (
	QObject,
	Signal,
	Slot,
	QCoreApplication,
	QTimer,
	QThread,
	QEventLoop,
	QMutex,
	QMutexLocker,
	QWaitCondition,
	)
from PySide2.QtWidgets import (
	QWidget,
	QDialog,
	QSizePolicy,
	QVBoxLayout,
	QHBoxLayout,
	QProgressBar,
	)

# 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 test")
db = doc.getDatabase(1)

# create evaluation options
# here we want to extract data for all steps of the last stage
all_stages = db.getStageIDs()
last_stage = all_stages[-1]
all_steps = db.getStepIDs(last_stage)
opt = MpcOdbVirtualResultEvaluationOptions()
opt.stage = last_stage

# create empty lists to store the maximum and minimum values for each node
env_max = []
env_min = []


class Worker(QObject):
	# signals
	finished = Signal()
	progress = Signal(int)
	# lengthy function
	@Slot()
	def run(self):
		try:
			# Nodal result
			NN = [6502] # list of nodes to process
			# calculate total number of nodes to process
			total_nodes = len(NN)
			# for each node
			for idx, node_id in enumerate(NN):
				# create empty list to store the acceleration values for this node
				acc_values = []
				# for each step
				for i in all_steps:
					opt.step = i
					# evaluate acceleration in x direction for this node
					field = db.getNodalResult("Acceleration", match=MpcOdb.Contains)
					acc_x = field[MpcOdbResultField.node(node_id), 0]
					# append the acceleration value to the list
					acc_values.append(acc_x)
				# calculate the maximum and minimum values for this node
				node_env_max = np.max(acc_values)
				node_env_min = np.min(acc_values)
				# append the maximum and minimum values to the respective lists
				env_max.append(node_env_max)
				env_min.append(node_env_min)
				# update progress bar
				self.progress.emit(int((idx + 1) / total_nodes * 100))

			# create virtual results
			max_vr = vr(env_max, 'Acceleration', ['Ax'], displayName='(ENVmax) Acceleration')
			min_vr = vr(env_min, 'Acceleration', ['Ax'], displayName='(ENVmin) Acceleration')

			# add virtual results to the database
			db.addVirtualResult(max_vr)
			db.addVirtualResult(min_vr)

		except Exception as ex:
			IO.write_cerr(str(ex))
		finally:
			self.finished.emit()
			
# create the thread, the worker, and move the worker to the thread
thread = QThread()
worker = Worker()
worker.moveToThread(thread)

# set up worker and thread connections
worker.finished.connect(thread.quit)
worker.finished.connect(worker.deleteLater)
thread.started.connect(worker.run)
thread.finished.connect(thread.deleteLater)

# create dialog
dialog = QDialog()
dialog.setWindowTitle("Processing nodes...")
#dialog.setWindowModality(2)
dialog.setFixedWidth(300)
dialog.setFixedHeight(100)
vbox = QVBoxLayout()
hbox = QHBoxLayout()
progress = QProgressBar()
hbox.addWidget(progress)
vbox.addLayout(hbox)
dialog.setLayout(vbox)

# create the thread, the worker,
# start the thread
thread.start()

# run dialog
dialog.exec_()

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

Re: Envelope of accelerations in OpenseesMP with API

Post by STKO Team » Fri Apr 21, 2023 7:45 am

There are some mistakes in your code.
But first I think you are misunderstanding how to use a virtual result.
If you use the MpcOdbVirtualResult.max (or min, or any other operator between virtual results), the output will be a virtual result that you can visualize in a plot. However, the operation is among the operands at the same time step.
I.e. MpcOdbVirtualResult.max(R1, R2) creates a result that gives you, for each time step, and for each node/element, the max between the operands R1 and R2.

Instead, I think you want the envelope in time. If this is the case, you have to evaluate a result, for each time step, create the field, and extract data. Then you can compare this data with the previous time step.

I think it's better if you explain in more details what result you are trying to obtain

gabrifio
Posts: 18
Joined: Wed Dec 16, 2020 12:35 pm

Re: Envelope of accelerations in OpenseesMP with API

Post by gabrifio » Fri Apr 21, 2023 11:35 am

Dear Massimo,

yes, you are correct. I am trying to plot the envelope of the maximum (and minimum) accelerations in time for each node. I would like to see a colour map like the Surface Color map but of course without time steps, only the maximum (and minimum)value for each node.

I have corrected the error and now the code seem to work until the point when I print "DONE". In total I have 6502 nodes but here I just put 30 nodes to speed up the calculation. Also, I have 8190 time steps but I am just plotting between 1555 and 1600 also to save time (it takes the strongest part of the accelerogram).

Here is the section of the code that works:

Code: Select all

from PyMpc import * 
from PyMpc import MpcOdbVirtualResult as vr
from time import sleep
import traceback
import numpy as np
from PySide2.QtCore import Qt, QObject, Signal, Slot, QThread, QEventLoop
from PySide2.QtWidgets import QDialog, QLabel, QProgressBar, QVBoxLayout, QApplication, QInputDialog
from numpy import genfromtxt
from numpy import argwhere
from numpy import sqrt

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 test")
db = doc.getDatabase(1)

# create evaluation options
# here we want to extract data for all steps of the last stage
all_stages = db.getStageIDs()
last_stage = all_stages[-1]
all_steps = db.getStepIDs(last_stage)
opt = MpcOdbVirtualResultEvaluationOptions()
opt.stage = last_stage
# for updating every n steps
num_steps = len(all_steps)

# initialization results
ACCELERATIONS_ENV_P = [];
# The worker class. It will run the lengthy function on a working thread
# and emit signals for calling the connected slots on the main thread
# Accelerations
acc = db.getNodalResult("Acceleration", match=MpcOdb.Contains)
		
class Worker(QObject):
	# signals
	finished = Signal()
	# lengthy function
	@Slot()
	def run(self):
		try:
		# forloop for all the nodes in the database
			for j in range(30): # Total number of nodes is 6502 but to save time I want to try with 30 nodes first.
				print(j+1)
				tmp = [];
				for i in range(1555, 1600, 1): #Total length = 8190, just taking the strongest part of signal.
					opt.step = all_steps[i];
					field = acc.evaluate(opt);
					row = MpcOdbResultField.node(j+1) # the node field object is the row identifier in the result field
					col = 0 # 0 = Ax, 1 = Ay, 2 = Az (the column 0-based index, depends on the result)
					AX = field[row, col];
					tmp.append(AX);
					AX = 0;
				ACCELERATIONS_ENV_P.append(max(tmp));         
		except Exception:
			IO.write_cerr(traceback.format_exc())
		finally:
			self.finished.emit()

# create an event loop that will block the script till the end
loop = QEventLoop()

# create the thread, the worker, and move the worker to the thread
thread = QThread()
worker = Worker()
worker.moveToThread(thread)

# set up worker and thread connections
worker.finished.connect(thread.quit)
worker.finished.connect(worker.deleteLater)
thread.started.connect(worker.run)
thread.finished.connect(thread.deleteLater)
thread.finished.connect(loop.quit)

# start the thread
thread.start()

# run event loop
loop.exec_()
print("DONE")
Now what I am not able to do is to load the results on the editor. I know it should be something like this:

Code: Select all

D1 = db.getNodalResult("Acceleration", match=MpcOdb.Contains)
D1 = D1 + ACCELERATIONS_ENV_P;
D1.setDisplayName("Envelope Max")
D1.setComponentLabels(["Envelope Max"])
db.addVirtualResult(D1)
But I get this error :
Cannot add the virtual results due to the following errors:
("Acceleration" + [0.0954871, 0.0484654, 0.0484654, 0.0954871, 0.0463059, 0.0463059, 0.0242182, 0.0242182, 0.0954288, 0.0951113, 0.094517, 0.0936013, 0.0923038, 0.0905579, 0.0883062, 0.0855173, 0.082202, 0.0784229, 0.0742939, 0.0699697, 0.0656259, 0.0614355, 0.0575467, 0.0540648, 0.0510415, 0.0486235, 0.0487237, 0.048782, 0.0488044, 0.0487954]) - Operands have different size (2 - 30)

I understand why I get this error but I don't know how to use db.addVirtualResult correctly to load the results in the STKO post-processor.

Thank you very much

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

Re: Envelope of accelerations in OpenseesMP with API

Post by STKO Team » Wed May 03, 2023 8:35 am

Unfortunately the time-envelope cannot be used as a virtual result.
In your case the result at time t, depends on all other previous time steps. Instead, a virtual result is an instantaneous result, i.e. an operation among multiple results at the same time step.

As I show in the webinar, you can extract the result for each node, and compute the envelope. but then you need to save those data as a file.

If you really want to visualize it in the post proc, you can do something like:
  • save the data in a text file
  • re-run a 1-step static analysis imposing the solution (disp and rot) to all DOFs, equal to the data you want to visualize.
  • plot the displacement of this fake analysis

gabrifio
Posts: 18
Joined: Wed Dec 16, 2020 12:35 pm

Re: Envelope of accelerations in OpenseesMP with API

Post by gabrifio » Thu May 04, 2023 9:30 am

Hi, thank you for the clarification. Actually, when I say envelope I would "only" need the maximum (positive) and minimum (negative) value reached in terms of acceleration for each node.

The problem is that (of course) the maxima are not simultaneous for all the nodes, so when you write:
Instead, a virtual result is an instantaneous result, i.e. an operation among multiple results at the same time step.
it seems to me that the most efficient way to do that is to extract the maximum and minimum acceleration for all the nodes and download it as a file, and after that plot the results outside the STKO post-processor. Of course the file should contain the coordinates of the node together with the acceleration values.

What do you think?

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

Re: Envelope of accelerations in OpenseesMP with API

Post by STKO Team » Wed May 17, 2023 3:42 pm

Yes now you can do it like this.

Post Reply