idaholab / mastodon

A MOOSE app for structural dynamics, seismic analysis, and risk assessment.
https://mooseframework.inl.gov/mastodon
GNU Lesser General Public License v2.1
40 stars 54 forks source link

FP Isolator : Post-yield stiffness terms (Lines 424-427) #333

Closed samyogshr closed 3 years ago

samyogshr commented 4 years ago

Could you please check calculation of stiffness terms (lines 424-427)? I have attached my calculation for K(1,1) term. Let me know whether this is an issue or not.

FP_stiffness.pdf

SaiSharathParsi1995 commented 4 years ago

@samyogshr and @cbolisetti , likely there is an issue. Please allow me sometime to dig into it.

samyogshr commented 4 years ago

@SaiSharathParsi1995 were you able to find out if it is a problem or not in the first place?

SaiSharathParsi1995 commented 4 years ago

@samyogshr, your calculations are correct and that is how they should be. The MASTODON element was adapted from the OpenSees FPBearingPTV element which somehow has the current equations with which I disagree as well. I have notified the developer of that element to look into it. Meanwhile we can try changing the stiffness equations in MASTODON to the correct ones and see how they influence the response. To my understanding, it shouldn't change anything significant. We can co-ordinate with @cbolisetti to work around this.

SaiSharathParsi1995 commented 4 years ago

Hello @samyogshr , I would recommend you to make the following changes to the code and see how the response changes for the MASTODON examples. Please verify if the new equations for stiffness terms are what you had derived.

Below are lines 415 to 428 in the original code: // Account for uplift Real k0Plastic = _k0;

  if (is_uplift)
    k0Plastic = kFactUplift * _k0;

  // Set tangent stiffness
  Real D = pow(qTrialNorm, 3);

  _Kb[_qp](1, 1) = k0Plastic * (1.0 - (qYield * _qTrial(1) * _qTrial(1) / D)) + k2y;
  _Kb[_qp](1, 2) = -qYield * k0Plastic * _qTrial(0) * _qTrial(1) / D;
  _Kb[_qp](2, 1) = _Kb[_qp](1, 2);
  _Kb[_qp](2, 2) = k0Plastic * (1.0 - (qYield * _qTrial(0) * _qTrial(0) / D)) + k2z;
}

Please change them to the following // Account for uplift

  if (is_uplift)
    k0y = kFactUplift * k0y;
    k0z = kFactUplift * k0z;

  // Set tangent stiffness
  Real D = pow(qTrialNorm, 3);

  _Kb[_qp](1, 1) = k0y * qYield * _qTrial(1) * _qTrial(1) / D + k2y;
  _Kb[_qp](1, 2) = -qYield * k0z * _qTrial(1) * _qTrial(0) / D;
  _Kb[_qp](2, 1) = -qYield * k0y * _qTrial(0) * _qTrial(1) / D;
  _Kb[_qp](2, 2) = k0z * qYield * _qTrial(0) * _qTrial(0) / D + k2z;
}   
samyogshr commented 4 years ago

Hello @SaiSharathParsi1995,

Could you please share how you derived Kb(1,1)? Another question - The term N*localrot will not show up, so how will that be accounted for?

Thank you Samyog

SaiSharathParsi1995 commented 4 years ago

@samyogshr, please see the screenshot attached. I just simplified the equation in your above pdf.

The term N*localrot has nothing to do with bearing deformation and hence doesn't appear in the stiffness term. It is simply additional force caused because of rotation of the concave plate as I mentioned in my email to you last week.

Screenshot_1

samyogshr commented 4 years ago

@SaiSharathParsi1995 Thank you and ko ~ koy ~ koz, so it shouldn't matter, I agree to this and am happy to accept this approximation. But, is this approximation necessary?

In relation to N*localrot term, isn't it true that : Fb(1,0) = Kb(1,1) basic_def(1,0) + Kb(1,2) basic_def(2,0)? Could you clarify this part?

Samyog

SaiSharathParsi1995 commented 4 years ago

@samyogshr , we are not making the approximation ko~koy~koz. They are different. The plasticity model is for the hysteretic part where the initial stiffness is koy and not ko.

ko is sum of koy and k2y. I don't see any approximation here.

In regard to you second question, in the expression for Kb (1,1) we already accounted for bi-directional effect. So, I am not sure if we have to add Kb(1,2)*basic_def(2,0). At this moment, I don't have a technical answer for this but will get back to you if I find one.

samyogshr commented 4 years ago

@SaiSharathParsi1995 In the off-diagonal terms Kb(1,2) and Kb(2,1), I am thinking that the code has "ko" term where it should be "koy" and "koz". So, if it is assumed that k2y << ko, ko ~ koy ~ koz, the expression in the code would be acceptable. But, this is not true for Kb(1,1) and Kb(2,2), which I didn't realize earlier.

Regarding the second question, I am thinking of the matrix form {Fb} = [Kb] {basic_def}, where off-diagonal terms in [kb] are non-zero once frictional resistance is overcome (plastic stage). Do you agree with this? If you do, there is no local rotation term in the right hand side of the matrix form, isn't it?

Samyog

SaiSharathParsi1995 commented 4 years ago

The off-diagonal terms should have koy and koz as seen in the modified equations I posted on this thread.

I agree with you. But I am pretty sure the formulation works since it had been in use for a long time now. Are the rotation terms being accounted when converting basic to local and global directions. This is just my assumption but we have to look into it.

cbolisetti commented 3 years ago

The formulation works for now. Will re-open when a bug is demonstrated.