Outputs for damping forces

Post Reply
DeltaGL
Posts: 9
Joined: Mon Oct 26, 2020 9:57 am

Outputs for damping forces

Post by DeltaGL » Thu Feb 22, 2024 5:53 pm

Hi,
I'm trying to obtain in a concentrated plasticity model the damping forces aiming at obtaining the damping energy by integration of damping forces and nodal displacements. I've used the model for damping "modalDamping" and the recorders "rayleighForces" and "rayleighForces!". With the first recorder, I obtained damping force equal to 0 in all nodes in each step and, with the second recorder, I obtained a line of data for each mode that I understand that they are the forces obtained for each mode and node in each step.

Could I use the recorder “rayleighForces” in "modalDamping" or it is only valid for rayleigh command?

Could someone explain me how I can obtain the total damping energy dissipated by the model along an earthquake? I mean the outputs that I would need to compute it.

I really appreciate any comment

Best regards and thanks in advance

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

Re: Outputs for damping forces

Post by STKO Team » Mon Feb 26, 2024 11:18 am

There is not an option to record damping coming from the modalDamping command as it works at the integrator level, and not at the element/node level as the rayleigh damping.

A simple hack you can use is to record the full right-hand-side (RHS) vector with printB at each step twice.
Before recording it the second time, you can turn off damping and reset it later on.
In this way the first printB will contain the damping forces while the second will not.

A simple python example to see if it fits your needs:

Code: Select all

from openseespy import opensees as ops
import numpy as np
import matplotlib.pyplot as plt

ops.wipe()
ops.model('basic', '-ndm', 3, '-ndf', 6)

ops.node(1, 0,0,0)
ops.node(2, 1,0,0, '-mass', 1,0,0,0,0,0)

ops.section('Elastic', 1, 1.0, 1, 1,1,1,1)
ops.geomTransf('Linear', 1,  0,0,1)
ops.beamIntegration('Lobatto', 1, 1, 3)
ops.element('dispBeamColumn', 1, 1,2, 1, 1)
ops.fix(1, 1,1,1, 1,1,1)
ops.fix(2, 0,1,1, 1,1,1)

ops.timeSeries('Path', 1,  '-dt', 0.01, '-values', *[0.0, 1.0, 0.0])
ops.pattern('Plain',1,1)
ops.load(2, 1,0,0, 0,0,0)

Nmodes = 1
w2 = ops.eigen('-fullGenLapack', Nmodes)
zeta = 0.2
damping = [zeta]*Nmodes
ops.modalDampingQ(*damping)

N = 1000
dt = 0.01
u = [0.0]*N
t = [0.0]*N
r = [0.0]*N
s = [0.0]*N
ctime = 0.0
ops.analysis('Transient')
ops.integrator('Newmark', 0.5, 0.25)
for i in range(N):
  ops.analyze(1, dt)
  ops.reactions()
  ctime = float(i+1)*dt
  t[i] = ctime
  u[i] = ops.nodeDisp(2)[0]
  r[i] = ops.nodeReaction(1)[0]
  node_rhs_with_damp = ops.printB('ret')[0]
  ops.modalDampingQ(*[0.0]*Nmodes)
  node_rhs_without_damp = ops.printB('ret')[0]
  s[i] = node_rhs_with_damp - node_rhs_without_damp
  ops.modalDampingQ(*damping)

plt.plot(t,r)
plt.plot(t,s,'--')
plt.grid(':')
plt.show()

DeltaGL
Posts: 9
Joined: Mon Oct 26, 2020 9:57 am

Re: Outputs for damping forces

Post by DeltaGL » Mon Feb 26, 2024 4:58 pm

Thank you very much! Good hack. Although it is not very good for large models because you need the fullGenLapack solver, slowing the process.

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

Re: Outputs for damping forces

Post by STKO Team » Tue Feb 27, 2024 2:11 pm

Actually I used the fullGenLapack simply because in that example I had only 1 dynamic DOF (1 free dof with mass) and asked for 1 eigenvalue.
The fast Arpack (default) solver can at most give you N-1 eigenvalues, where N is the number of dynamic DOFs.

In this example I use the printB which works also with sparse systems. What instead works only with dense solvers is the printA (which is not used here).

See this example in which I have 3 dynamic DOFs and ask for only 1 mode with Arpack and UmfPack (sparse solvers).

Code: Select all

from openseespy import opensees as ops
import numpy as np
import matplotlib.pyplot as plt

ops.wipe()
ops.model('basic', '-ndm', 3, '-ndf', 6)

ops.node(1, 0,0,0)
ops.node(2, 1,0,0, '-mass', 1,1,1,0,0,0)

ops.section('Elastic', 1, 1.0, 1, 1,1,1,1)
ops.geomTransf('Linear', 1,  0,0,1)
ops.beamIntegration('Lobatto', 1, 1, 3)
ops.element('dispBeamColumn', 1, 1,2, 1, 1)
ops.fix(1, 1,1,1, 1,1,1)
ops.fix(2, 0,0,0, 1,1,1)

ops.timeSeries('Path', 1,  '-dt', 0.01, '-values', *[0.0, 1.0, 0.0])
ops.pattern('Plain',1,1)
ops.load(2, 1,0,0, 0,0,0)

Nmodes = 1
w2 = ops.eigen(Nmodes)
zeta = 0.2
damping = [zeta]*Nmodes
ops.modalDampingQ(*damping)

N = 1000
dt = 0.01
u = [0.0]*N
t = [0.0]*N
r = [0.0]*N
s = [0.0]*N
ctime = 0.0
ops.analysis('Transient')
ops.system('UmfPack')
ops.integrator('Newmark', 0.5, 0.25)
for i in range(N):
  ops.analyze(1, dt)
  ops.reactions()
  ctime = float(i+1)*dt
  t[i] = ctime
  u[i] = ops.nodeDisp(2)[0]
  r[i] = ops.nodeReaction(1)[0]
  node_rhs_with_damp = ops.printB('ret')[0]
  ops.modalDampingQ(*[0.0]*Nmodes)
  node_rhs_without_damp = ops.printB('ret')[0]
  s[i] = node_rhs_with_damp - node_rhs_without_damp
  ops.modalDampingQ(*damping)

plt.plot(t,r)
plt.plot(t,s,'--')
plt.grid(':')
plt.show()



Post Reply