Parameter optimization for biphasic confined compression stress-relaxation

Collapse
X
 
  • Filter
  • Time
  • Show
Clear All
new posts
  • ateshian
    Developer
    • Dec 2007
    • 1821

    Parameter optimization for biphasic confined compression stress-relaxation

    The two attached files provide an example of how to perform a parameter optimization (curve-fitting) for confined compression stress-relaxation.

    CCsrlxOpt.feb is the input file for performing the parameter optimization. Execute from the prompt:
    Code:
    febio -s CCsrlxOpt.feb
    CCsrlxFit.feb is the finite element model for biphasic confined compression stress-relaxation.

    Model:
    Since confined compression is a one-dimensional analysis, only the height (thickness) of the material matters. Therefore a cuboidal model is used here (with a height of 1 mm).

    A rigid body is used to represent the indenter. The indenter is connected via a rigid interface to the top surface (z=1) of the cuboid. The fluid pressure is set to zero on that surface, because the indenter is idealized as free-draining. The indenter is given a prescribed displacement (ramp-and-hold).

    Optimization:
    The optimization is performed on the reaction force Fz of the indenter. In this example, the "experimental" data appearing in CCsrlxOpt.feb was actually generated in FEBio, therefore the parameter optimization should converge to exact values of E=1 MPa and perm=0.001 mm^4/N.s.

    The initial guess for E is taken to be 0.5 MPa and a representative scale for this parameter is 1. The initial guess for perm is 5e-4 mm^4/N.s and its representative scale is 1e-3. With this choice of initial guesses, the problem converges nicely.

    Note: For some values of the initial guess for perm (e.g., 0.002), the optimization algorithm may fail in this example, because it might produce a negative value for the permeability during the iterative optimization process. (As of the writing of this example, parameter optimization in FEBio does not yet implement bounds on the parameters.)

    Gerard
  • mjf2152
    Member
    • Jul 2011
    • 45

    #2
    Does this example require a certain version of FEBio? I get an error message when I try to run the optimization:
    FATAL ERROR: the variable Biphasic.solid.E is not recognized.
    The downloads don't seem to be working at the moment so I cannot upgrade, but I am currently using version 1.4.1.2762.

    After trying a few different things, it seems like my FEBio does not correctly descend the XML tree to the sub-property level. By that I mean that I can get it to recognize material.property, but not material.solid.property.

    Comment

    • ateshian
      Developer
      • Dec 2007
      • 1821

      #3
      Please use the latest FEBio development version (09-Jan-2012 as of this posting).
      Last edited by ateshian; 01-09-2012, 07:39 PM.

      Comment

      • mjf2152
        Member
        • Jul 2011
        • 45

        #4
        Thank you! It is working for me now.

        Comment

        • ntjacobs
          Member
          • Nov 2009
          • 79

          #5
          Hi Gerard,

          Thanks for posting this example. I've been trying to run the optimization and I've not been able to replicate your results of of E=1 MPa and perm=0.001 mm^4/N.s. I'm getting E=0.5 and perm = 0.0005. I was hoping you might have an idea of where I've gone wrong. I didn't modify either of the files you posted, and I ran the optimization from the command prompt as febio -s CCsrlxopt.feb. I've attached the .log and .dat files the run produced in the zipped folder NTJ_Optimization.zip.

          Thanks for your help,

          Nathan
          Attached Files

          Comment

          • ntjacobs
            Member
            • Nov 2009
            • 79

            #6
            I forgot to mention, I've been using FEBio for WinXP version 1.5.1.3739

            Comment

            • ateshian
              Developer
              • Dec 2007
              • 1821

              #7
              Hi Nathan,

              I've retried these files on my Intel Mac (64 bit) and I get the expected convergence (E=1 MPa and k=0.001 mm^4/N.s). The parameter optimization also works on a Windows 7 64 bit machine. It would seem to me that perhaps the problem might lie with the WinXP version? If that's the case I would need to defer to some of my colleagues on the FEBio team, because it would not be a coding issue as much as an architecture issue. I'll ask Dave to look into it.

              Best,

              Gerard

              Comment

              • ateshian
                Developer
                • Dec 2007
                • 1821

                #8
                Hi Nathan,

                Right after I posted my response below, one of my students informed me he was able to run this optimization problem successfully on a Windows XP machine in our lab. (He downloaded the latest version of FEBio). So I am not sure that the problem is related to the architecture. Did you compile FEBio yourself or are you using the executable from the software download website? Would you be able to try this analysis on another computer?

                Best,

                Gerard

                Comment

                • Farrell1
                  Junior Member
                  • Oct 2011
                  • 25

                  #9
                  I have run into issues when fitting a non-linear biphasic model to multiple strain degenerate nucleus pulposus confined compression data.

                  I have manually obtained E and beta by fitting the non-linear biphasic equation to equilibrium stress data thus reducing the parameter space when using FEBio to fit permeability and M, however, I can't seem to obtain a good transient fit. I have also tried simplifying the model to constant/isotropic permeability as well as leaving all parameters open to fit through FEBio, neither of which yielded improvement.

                  I have attached both of my files for comment incase I am doing something wrong.

                  Thanks in advance,

                  Mark
                  Attached Files

                  Comment

                  • ateshian
                    Developer
                    • Dec 2007
                    • 1821

                    #10
                    Hi Mark,

                    Looking at your data, it does indeed appear that the behavior observed in your experiments cannot be predicted by a biphasic analysis with a Holmes-Mow solid matrix and strain-dependent permeability. There are a couple of potential ways to address this issue.

                    I notice that the load response to the first ramp behaves substantially differently from that of the second and third ramps: The first ramp shows a very small peak-to-equilibrium ratio compared to the other two; it also equilibrates to a steady-state value much faster than the other two. If this is a true behavior of the tissue in confined compression, then I doubt that the model will be able to capture this effect no matter what choice of material properties you select.

                    On the other hand, it is possible that the tissue sample was not fully confined in this first ramp, so that its behavior isn't fully consistent with predictions from a confined compression analysis. (In my experience with cartilage, even if you size the sample to match the chamber diameter exactly, and even if the sample swells a little bit inside the chamber, you still need to apply a substantial tare strain, e.g., up to 15%, before the response becomes consistent with theoretical confined compression.) You can verify this possibility by choosing to treat the response to the first ramp as a tare condition and fitting only the second and third ramps (with the load and displacement suitably shifted to start at the origin). If the fit improves considerably, this becomes a plausible explanation.

                    Another possibility is to increase the ramp time in the experiment. Currently, your results for the 2nd and 3rd ramps show that relaxation occurs much more slowly than the ramp time (you ramp up the displacement in ~50 s but it takes ~2000 s for the load to relax). It is possible that this protocol excites the intrinsic viscoelasticity of the solid matrix (which may or may not be linear viscoelasticity). So it might be helpful to try ramping up more slowly (e.g., in 500 s instead of 50 s) and seeing whether that produces better fits.

                    Let me know if these suggestions help.

                    Best,

                    Gerard

                    Comment

                    • Farrell1
                      Junior Member
                      • Oct 2011
                      • 25

                      #11
                      Gerard,

                      Thank you for the informative reply. Your thoughts are similar to my conclusions but I wanted to check that I was running the anlysis correctly. I will try to fit to the 2nd and 3rd ramps only as I believe the first ramp is revealing an experimental artefact as you suggest - the degenerate tissue was significantly depleted of PG content thus negligible swelling occurred during the initial hold phase before the onset of deformation so the tissue may not have occupied the full chamber volume.

                      Kind regard,

                      Mark

                      Comment

                      • ateshian
                        Developer
                        • Dec 2007
                        • 1821

                        #12
                        Hi Mark,

                        Yes, your optimization analysis is set up properly. However, you can improve the efficiency of the computations in two ways:

                        - For a stress-relaxation analysis with a ramp time of 50 s and relaxation time of ~2000 s, I don't recommend using constant time increments. You should create a must-point curve that controls the maximum time increment (dtmax), e.g., as follows: 5 s max increments from 0 to 50 s, 5 s max increments from 50 to 100 s, then 100 s max increments from 100 s to the end of the relaxation. Repeat for each stress-relaxation step. Make sure your time points correspond exactly to the beginning and end of each ramp.

                        - I have found that optimizations involving the hydraulic permeability converge faster if the parameter <f_diff_scale> is set to 0.01 (instead of the default of 0.001) in the optimization input file. That's because very small changes in the permeability don't affect the solution as much.

                        These are the settings I used with your files to increase computational efficiency.

                        Best,

                        Gerard

                        Comment

                        • hamed71200
                          Junior Member
                          • Apr 2013
                          • 24

                          #13
                          Dr. Ateshian

                          I am working on a biphasic model and I want to use a "Transversely isotropic viscoelastic" solid phase.
                          I know how to run the parameter optimization for the "Isotropic viscoelastic" solid phase.
                          My problem is that I could not find "Transversely isotropic viscoelastic" option.
                          I think I should choose orthotropic viscoelastic and make two of elastic moduli, shear moduli and Poisson's ratios equal to each other. Is this right?

                          I would be grateful if you could help me how to define the transverse isotropic viscoleasticity. if what I said above it's right, how should put the moduli equal in the optimization file?

                          Thanks,
                          Hamed

                          Comment

                          • ateshian
                            Developer
                            • Dec 2007
                            • 1821

                            #14
                            Hi Hamed,

                            I think I should choose orthotropic viscoelastic and make two of elastic moduli, shear moduli and Poisson's ratios equal to each other. Is this right?
                            Yes, that's correct.

                            if what I said above it's right, how should put the moduli equal in the optimization file?
                            Unfortunately it is not yet possible to constrain multiple properties to have the same value during an optimization. You can submit a feature request to add a "transversely isotropic" material in FEBio, which is easier to do than implementing constrained optimization.

                            Best,

                            Gerard

                            Comment

                            • hamed71200
                              Junior Member
                              • Apr 2013
                              • 24

                              #15
                              Thanks Dr. Ateshian

                              Originally posted by ateshian View Post
                              You can submit a feature request to add a "transversely isotropic" material in FEBio
                              I did not understand this part. How can I submit a feature request to add "transversely isotropic" ?

                              Regards,

                              Hamed

                              Comment

                              Working...
                              X
                              😀
                              😂
                              🥰
                              😘
                              🤢
                              😎
                              😞
                              😡
                              👍
                              👎