.. DO NOT EDIT.
.. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY.
.. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE:
.. "content/examples/10-pgi/plot_inv_1_PGI_Linear_1D_joint_WithRelationships.py"
.. LINE NUMBERS ARE GIVEN BELOW.

.. only:: html

    .. note::
        :class: sphx-glr-download-link-note

        :ref:`Go to the end <sphx_glr_download_content_examples_10-pgi_plot_inv_1_PGI_Linear_1D_joint_WithRelationships.py>`
        to download the full example code

.. rst-class:: sphx-glr-example-title

.. _sphx_glr_content_examples_10-pgi_plot_inv_1_PGI_Linear_1D_joint_WithRelationships.py:


Petrophysically guided inversion: Joint linear example with nonlinear relationships
===================================================================================

We do a comparison between the classic least-squares inversion
and our formulation of a petrophysically guided inversion.
We explore it through coupling two linear problems whose respective physical
properties are linked by polynomial relationships that change between rock units.

.. GENERATED FROM PYTHON SOURCE LINES 11-432



.. image-sg:: /content/examples/10-pgi/images/sphx_glr_plot_inv_1_PGI_Linear_1D_joint_WithRelationships_001.png
   :alt: Problem 1, Problem 2, Petrophysical Distribution, Problem 1, Problem 2, Petrophysical Distribution, Problem 1, Problem 2, Petro Distribution
   :srcset: /content/examples/10-pgi/images/sphx_glr_plot_inv_1_PGI_Linear_1D_joint_WithRelationships_001.png
   :class: sphx-glr-single-img


.. rst-class:: sphx-glr-script-out

 .. code-block:: none

    SimPEG.InvProblem will set Regularization.reference_model to m0.
    SimPEG.InvProblem will set Regularization.reference_model to m0.
    SimPEG.InvProblem will set Regularization.reference_model to m0.

                        SimPEG.InvProblem is setting bfgsH0 to the inverse of the eval2Deriv.
                        ***Done using the default solver Pardiso and no solver_opts.***
                    
    Alpha scales: [3.4879629923623376, 0.0, 3.4844686613163604e-06, 0.0]
    Calculating the scaling parameter.
    Scale Multipliers:  [0.09605361 0.90394639]
    <class 'SimPEG.regularization.pgi.PGIsmallness'>
    Initial data misfit scales:  [0.09605361 0.90394639]
    model has any nan: 0
    =============================== Projected GNCG ===============================
      #     beta     phi_d     phi_m       f      |proj(x-g)-x|  LS    Comment   
    -----------------------------------------------------------------------------
    x0 has any nan: 0
       0  1.86e+01  1.50e+05  0.00e+00  1.50e+05    1.41e+02      0              
    geophys. misfits: 498.1 (target 15.0 [False]); 35.7 (target 15.0 [False]) | smallness misfit: 1459.2 (target: 100.0 [False])
    Beta cooling evaluation: progress: [498.1  35.7]; minimum progress targets: [120000. 120000.]
       1  1.86e+01  8.01e+01  2.06e+01  4.65e+02    7.96e+01      0              
    geophys. misfits: 235.0 (target 15.0 [False]); 10.9 (target 15.0 [True]) | smallness misfit: 644.7 (target: 100.0 [False])
    Beta cooling evaluation: progress: [235.   10.9]; minimum progress targets: [398.5  28.6]
    Updating scaling for data misfits by  1.374269734362692
    New scales: [0.12742273 0.87257727]
       2  1.86e+01  3.95e+01  2.00e+01  4.11e+02    5.66e+01      0   Skip BFGS  
    geophys. misfits: 156.2 (target 15.0 [False]); 10.9 (target 15.0 [True]) | smallness misfit: 609.7 (target: 100.0 [False])
    Beta cooling evaluation: progress: [156.2  10.9]; minimum progress targets: [188.  15.]
    Updating scaling for data misfits by  1.381664440666124
    New scales: [0.16789048 0.83210952]
       3  1.86e+01  3.53e+01  2.04e+01  4.15e+02    5.64e+01      0   Skip BFGS  
    geophys. misfits: 105.3 (target 15.0 [False]); 10.9 (target 15.0 [True]) | smallness misfit: 586.1 (target: 100.0 [False])
    Beta cooling evaluation: progress: [105.3  10.9]; minimum progress targets: [125.  15.]
    Updating scaling for data misfits by  1.3811013876658975
    New scales: [0.21792991 0.78207009]
       4  1.86e+01  3.15e+01  2.08e+01  4.19e+02    5.71e+01      0              
    geophys. misfits: 73.0 (target 15.0 [False]); 11.0 (target 15.0 [True]) | smallness misfit: 568.0 (target: 100.0 [False])
    Beta cooling evaluation: progress: [73. 11.]; minimum progress targets: [84.3 15. ]
    Updating scaling for data misfits by  1.368561549136323
    New scales: [0.27607591 0.72392409]
       5  1.86e+01  2.81e+01  2.11e+01  4.22e+02    6.02e+01      0   Skip BFGS  
    geophys. misfits: 52.7 (target 15.0 [False]); 11.1 (target 15.0 [True]) | smallness misfit: 553.5 (target: 100.0 [False])
    Beta cooling evaluation: progress: [52.7 11.1]; minimum progress targets: [58.4 15. ]
    Updating scaling for data misfits by  1.34723481121236
    New scales: [0.33940285 0.66059715]
       6  1.86e+01  2.53e+01  2.14e+01  4.24e+02    5.20e+01      0   Skip BFGS  
    geophys. misfits: 40.2 (target 15.0 [False]); 11.4 (target 15.0 [True]) | smallness misfit: 541.1 (target: 100.0 [False])
    Beta cooling evaluation: progress: [40.2 11.4]; minimum progress targets: [42.2 15. ]
    Updating scaling for data misfits by  1.3174153484248725
    New scales: [0.40364876 0.59635124]
       7  1.86e+01  2.30e+01  2.16e+01  4.25e+02    4.71e+01      0   Skip BFGS  
    geophys. misfits: 32.4 (target 15.0 [False]); 11.7 (target 15.0 [True]) | smallness misfit: 530.2 (target: 100.0 [False])
    Beta cooling evaluation: progress: [32.4 11.7]; minimum progress targets: [32.2 15. ]
    Decreasing beta to counter data misfit decrase plateau.
    Updating scaling for data misfits by  1.27930152260321
    New scales: [0.46406942 0.53593058]
       8  9.31e+00  2.13e+01  2.17e+01  2.24e+02    7.29e+01      0   Skip BFGS  
    geophys. misfits: 14.1 (target 15.0 [True]); 9.7 (target 15.0 [True]) | smallness misfit: 556.7 (target: 100.0 [False])
    Beta cooling evaluation: progress: [14.1  9.7]; minimum progress targets: [25.9 15. ]
    Warming alpha_pgi to favor clustering:  1.3060218468359308
       9  9.31e+00  1.17e+01  2.28e+01  2.24e+02    2.11e+01      0              
    geophys. misfits: 13.9 (target 15.0 [True]); 10.3 (target 15.0 [True]) | smallness misfit: 515.5 (target: 100.0 [False])
    Beta cooling evaluation: progress: [13.9 10.3]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  1.6509719957323006
      10  9.31e+00  1.20e+01  2.31e+01  2.27e+02    1.78e+01      0              
    geophys. misfits: 13.8 (target 15.0 [True]); 10.9 (target 15.0 [True]) | smallness misfit: 481.9 (target: 100.0 [False])
    Beta cooling evaluation: progress: [13.8 10.9]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  2.0352960333864827
      11  9.31e+00  1.22e+01  2.35e+01  2.31e+02    3.67e+01      0   Skip BFGS  
    geophys. misfits: 13.7 (target 15.0 [True]); 11.5 (target 15.0 [True]) | smallness misfit: 451.6 (target: 100.0 [False])
    Beta cooling evaluation: progress: [13.7 11.5]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  2.440521801517239
      12  9.31e+00  1.25e+01  2.38e+01  2.34e+02    1.95e+01      0              
    geophys. misfits: 13.6 (target 15.0 [True]); 12.2 (target 15.0 [True]) | smallness misfit: 424.9 (target: 100.0 [False])
    Beta cooling evaluation: progress: [13.6 12.2]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  2.848206268636033
      13  9.31e+00  1.28e+01  2.42e+01  2.38e+02    1.93e+01      0   Skip BFGS  
    geophys. misfits: 13.6 (target 15.0 [True]); 12.8 (target 15.0 [True]) | smallness misfit: 402.0 (target: 100.0 [False])
    Beta cooling evaluation: progress: [13.6 12.8]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  3.2390019898577744
      14  9.31e+00  1.32e+01  2.45e+01  2.41e+02    2.84e+01      0   Skip BFGS  
    geophys. misfits: 13.6 (target 15.0 [True]); 13.5 (target 15.0 [True]) | smallness misfit: 382.8 (target: 100.0 [False])
    Beta cooling evaluation: progress: [13.6 13.5]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  3.5973567509010813
      15  9.31e+00  1.35e+01  2.48e+01  2.44e+02    4.30e+01      0   Skip BFGS  
    geophys. misfits: 13.6 (target 15.0 [True]); 14.0 (target 15.0 [True]) | smallness misfit: 367.0 (target: 100.0 [False])
    Beta cooling evaluation: progress: [13.6 14. ]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  3.9132521705024024
      16  9.31e+00  1.38e+01  2.50e+01  2.46e+02    2.11e+01      0   Skip BFGS  
    geophys. misfits: 13.6 (target 15.0 [True]); 14.5 (target 15.0 [True]) | smallness misfit: 354.5 (target: 100.0 [False])
    Beta cooling evaluation: progress: [13.6 14.5]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  4.182578893798094
      17  9.31e+00  1.41e+01  2.52e+01  2.48e+02    1.22e+01      0   Skip BFGS  
    geophys. misfits: 13.6 (target 15.0 [True]); 14.9 (target 15.0 [True]) | smallness misfit: 344.5 (target: 100.0 [False])
    Beta cooling evaluation: progress: [13.6 14.9]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  4.404910734903503
      18  9.31e+00  1.43e+01  2.53e+01  2.50e+02    1.02e+01      0   Skip BFGS  
    geophys. misfits: 13.6 (target 15.0 [True]); 15.3 (target 15.0 [False]) | smallness misfit: 336.8 (target: 100.0 [False])
    Beta cooling evaluation: progress: [13.6 15.3]; minimum progress targets: [15. 15.]
    Decreasing beta to counter data misfit increase.
    Updating scaling for data misfits by  1.1003745588507636
    New scales: [0.44037973 0.55962027]
      19  4.65e+00  1.46e+01  2.53e+01  1.32e+02    7.38e+01      0   Skip BFGS  
    geophys. misfits: 9.9 (target 15.0 [True]); 10.9 (target 15.0 [True]) | smallness misfit: 369.3 (target: 100.0 [False])
    Beta cooling evaluation: progress: [ 9.9 10.9]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  6.387706037636578
      20  4.65e+00  1.04e+01  2.75e+01  1.38e+02    3.69e+01      0              
    geophys. misfits: 10.0 (target 15.0 [True]); 11.7 (target 15.0 [True]) | smallness misfit: 313.3 (target: 100.0 [False])
    Beta cooling evaluation: progress: [10.  11.7]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  8.867527960511257
      21  4.65e+00  1.10e+01  2.90e+01  1.46e+02    6.39e+01      0              
    geophys. misfits: 10.3 (target 15.0 [True]); 12.8 (target 15.0 [True]) | smallness misfit: 267.8 (target: 100.0 [False])
    Beta cooling evaluation: progress: [10.3 12.8]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  11.61834393430031
      22  4.65e+00  1.17e+01  3.03e+01  1.53e+02    8.26e+01      0              
    geophys. misfits: 10.5 (target 15.0 [True]); 13.1 (target 15.0 [True]) | smallness misfit: 241.4 (target: 100.0 [False])
    Beta cooling evaluation: progress: [10.5 13.1]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  14.926843025096527
      23  4.65e+00  1.20e+01  3.20e+01  1.61e+02    7.42e+01      0              
    geophys. misfits: 11.2 (target 15.0 [True]); 15.1 (target 15.0 [False]) | smallness misfit: 199.9 (target: 100.0 [False])
    Beta cooling evaluation: progress: [11.2 15.1]; minimum progress targets: [15. 15.]
    Decreasing beta to counter data misfit increase.
    Updating scaling for data misfits by  1.3442902520887736
    New scales: [0.36923798 0.63076202]
      24  2.33e+00  1.37e+01  3.15e+01  8.70e+01    7.61e+01      0              
    geophys. misfits: 9.9 (target 15.0 [True]); 11.0 (target 15.0 [True]) | smallness misfit: 224.8 (target: 100.0 [False])
    Beta cooling evaluation: progress: [ 9.9 11. ]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  21.50041230040881
      25  2.33e+00  1.06e+01  3.57e+01  9.36e+01    8.13e+01      0              
    geophys. misfits: 10.3 (target 15.0 [True]); 11.2 (target 15.0 [True]) | smallness misfit: 185.3 (target: 100.0 [False])
    Beta cooling evaluation: progress: [10.3 11.2]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  30.12426037525881
      26  2.33e+00  1.08e+01  3.88e+01  1.01e+02    9.53e+01      0              
    geophys. misfits: 11.6 (target 15.0 [True]); 10.6 (target 15.0 [True]) | smallness misfit: 170.0 (target: 100.0 [False])
    Beta cooling evaluation: progress: [11.6 10.6]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  40.864544419405185
      27  2.33e+00  1.09e+01  4.25e+01  1.10e+02    9.70e+01      0              
    geophys. misfits: 11.3 (target 15.0 [True]); 10.1 (target 15.0 [True]) | smallness misfit: 144.0 (target: 100.0 [False])
    Beta cooling evaluation: progress: [11.3 10.1]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  57.7058247954592
      28  2.33e+00  1.05e+01  4.73e+01  1.21e+02    1.04e+02      0              
    geophys. misfits: 12.6 (target 15.0 [True]); 11.4 (target 15.0 [True]) | smallness misfit: 116.0 (target: 100.0 [False])
    Beta cooling evaluation: progress: [12.6 11.4]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  72.30312813814554
      29  2.33e+00  1.18e+01  4.95e+01  1.27e+02    1.04e+02      0              
    geophys. misfits: 18.3 (target 15.0 [False]); 15.5 (target 15.0 [False]) | smallness misfit: 98.6 (target: 100.0 [True])
    Beta cooling evaluation: progress: [18.3 15.5]; minimum progress targets: [15. 15.]
    Decreasing beta to counter data misfit increase.
      30  1.16e+00  1.66e+01  4.66e+01  7.07e+01    9.94e+01      0              
    geophys. misfits: 12.3 (target 15.0 [True]); 10.9 (target 15.0 [True]) | smallness misfit: 115.4 (target: 100.0 [False])
    Beta cooling evaluation: progress: [12.3 10.9]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  93.78457997050192
      31  1.16e+00  1.14e+01  5.60e+01  7.66e+01    9.25e+01      0              
    geophys. misfits: 10.9 (target 15.0 [True]); 10.2 (target 15.0 [True]) | smallness misfit: 109.1 (target: 100.0 [False])
    Beta cooling evaluation: progress: [10.9 10.2]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  133.62216327456488
      32  1.16e+00  1.04e+01  6.64e+01  8.77e+01    1.03e+02      0              
    geophys. misfits: 11.8 (target 15.0 [True]); 12.0 (target 15.0 [True]) | smallness misfit: 91.0 (target: 100.0 [True])
    All targets have been reached
    Beta cooling evaluation: progress: [11.8 12. ]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  168.4889550423769
    ------------------------- STOP! -------------------------
    1 : |fc-fOld| = 0.0000e+00 <= tolF*(1+|f0|) = 1.5000e+04
    0 : |xc-x_last| = 4.3858e-01 <= tolX*(1+|x0|) = 1.0000e-06
    0 : |proj(x-g)-x|    = 1.0274e+02 <= tolG          = 1.0000e-01
    0 : |proj(x-g)-x|    = 1.0274e+02 <= 1e3*eps       = 1.0000e-02
    0 : maxIter   =      50    <= iter          =     33
    ------------------------- DONE! -------------------------
    SimPEG.InvProblem will set Regularization.reference_model to m0.
    SimPEG.InvProblem will set Regularization.reference_model to m0.
    SimPEG.InvProblem will set Regularization.reference_model to m0.

                        SimPEG.InvProblem is setting bfgsH0 to the inverse of the eval2Deriv.
                        ***Done using the default solver Pardiso and no solver_opts.***
                    
    Alpha scales: [0.0003484547914829826, 0.0, 3.564944499171368e-06, 0.0]
    Calculating the scaling parameter.
    Scale Multipliers:  [0.09605361 0.90394639]
    <class 'SimPEG.regularization.pgi.PGIsmallness'>
    Initial data misfit scales:  [0.09605361 0.90394639]
    model has any nan: 0
    =============================== Projected GNCG ===============================
      #     beta     phi_d     phi_m       f      |proj(x-g)-x|  LS    Comment   
    -----------------------------------------------------------------------------
    x0 has any nan: 0
       0  1.97e+03  1.50e+05  0.00e+00  1.50e+05    1.41e+02      0              
    geophys. misfits: 44567.6 (target 15.0 [False]); 31701.9 (target 15.0 [False]) | smallness misfit: 145.4 (target: 100.0 [False])
    Beta cooling evaluation: progress: [44567.6 31701.9]; minimum progress targets: [120000. 120000.]
       1  1.97e+03  3.29e+04  3.41e-01  3.36e+04    8.98e+01      0              
    geophys. misfits: 329.7 (target 15.0 [False]); 11.9 (target 15.0 [True]) | smallness misfit: 53.1 (target: 100.0 [True])
    Beta cooling evaluation: progress: [329.7  11.9]; minimum progress targets: [35654.1 25361.5]
    Updating scaling for data misfits by  1.2588248712215981
    New scales: [0.11798153 0.88201847]
       2  1.97e+03  4.94e+01  1.53e-01  3.52e+02    8.08e+01      0   Skip BFGS  
    geophys. misfits: 50.3 (target 15.0 [False]); 12.8 (target 15.0 [True]) | smallness misfit: 23.3 (target: 100.0 [True])
    Beta cooling evaluation: progress: [50.3 12.8]; minimum progress targets: [263.7  15. ]
    Updating scaling for data misfits by  1.1684921475591683
    New scales: [0.13517339 0.86482661]
       3  1.97e+03  1.79e+01  6.52e-02  1.47e+02    7.91e+01      0   Skip BFGS  
    geophys. misfits: 42.5 (target 15.0 [False]); 10.5 (target 15.0 [True]) | smallness misfit: 23.0 (target: 100.0 [True])
    Beta cooling evaluation: progress: [42.5 10.5]; minimum progress targets: [40.3 15. ]
    Decreasing beta to counter data misfit decrase plateau.
    Updating scaling for data misfits by  1.4221529080495827
    New scales: [0.18185962 0.81814038]
       4  9.87e+02  1.64e+01  6.31e-02  7.87e+01    8.36e+01      0              
    geophys. misfits: 15.4 (target 15.0 [False]); 9.0 (target 15.0 [True]) | smallness misfit: 24.6 (target: 100.0 [True])
    Beta cooling evaluation: progress: [15.4  9. ]; minimum progress targets: [34. 15.]
    Updating scaling for data misfits by  1.663521297584011
    New scales: [0.26995276 0.73004724]
       5  9.87e+02  1.07e+01  6.68e-02  7.67e+01    6.49e+01      0              
    geophys. misfits: 12.0 (target 15.0 [True]); 9.4 (target 15.0 [True]) | smallness misfit: 24.8 (target: 100.0 [True])
    All targets have been reached
    Beta cooling evaluation: progress: [12.   9.4]; minimum progress targets: [15. 15.]
    Warming alpha_pgi to favor clustering:  1.4238182762757154
    ------------------------- STOP! -------------------------
    1 : |fc-fOld| = 0.0000e+00 <= tolF*(1+|f0|) = 1.5000e+04
    0 : |xc-x_last| = 5.7836e-02 <= tolX*(1+|x0|) = 1.0000e-06
    0 : |proj(x-g)-x|    = 6.4921e+01 <= tolG          = 1.0000e-01
    0 : |proj(x-g)-x|    = 6.4921e+01 <= 1e3*eps       = 1.0000e-02
    0 : maxIter   =      50    <= iter          =      6
    ------------------------- DONE! -------------------------
    SimPEG.InvProblem will set Regularization.reference_model to m0.
    SimPEG.InvProblem will set Regularization.reference_model to m0.

                        SimPEG.InvProblem is setting bfgsH0 to the inverse of the eval2Deriv.
                        ***Done using the default solver Pardiso and no solver_opts.***
                    
    Alpha scales: [3.511758442591592e-05, 0.0, 3.5141198769325145e-05, 0.0]
    Calculating the scaling parameter.
    Scale Multipliers:  [0.09605361 0.90394639]
    /home/vsts/work/1/s/SimPEG/directives/directives.py:334: UserWarning:

    There is no PGI regularization. Smallness target is turned off (TriggerSmall flag)

    Initial data misfit scales:  [0.09605361 0.90394639]
    model has any nan: 0
    =============================== Projected GNCG ===============================
      #     beta     phi_d     phi_m       f      |proj(x-g)-x|  LS    Comment   
    -----------------------------------------------------------------------------
    x0 has any nan: 0
       0  1.02e+06  1.50e+05  0.00e+00  1.50e+05    1.41e+02      0              
    geophys. misfits: 27709.3 (target 15.0 [False]); 17519.4 (target 15.0 [False])
       1  2.05e+05  1.85e+04  2.16e-02  2.29e+04    1.37e+02      0              
    geophys. misfits: 4050.6 (target 15.0 [False]); 1865.0 (target 15.0 [False])
       2  4.10e+04  2.07e+03  5.34e-02  4.26e+03    1.24e+02      0   Skip BFGS  
    geophys. misfits: 268.2 (target 15.0 [False]); 118.9 (target 15.0 [False])
       3  8.20e+03  1.33e+02  7.01e-02  7.08e+02    9.19e+01      0   Skip BFGS  
    geophys. misfits: 20.6 (target 15.0 [False]); 15.5 (target 15.0 [False])
       4  1.64e+03  1.60e+01  7.51e-02  1.39e+02    6.28e+01      0   Skip BFGS  
    geophys. misfits: 9.2 (target 15.0 [True]); 8.1 (target 15.0 [True])
    All targets have been reached
    ------------------------- STOP! -------------------------
    1 : |fc-fOld| = 0.0000e+00 <= tolF*(1+|f0|) = 1.5000e+04
    0 : |xc-x_last| = 4.3431e-01 <= tolX*(1+|x0|) = 1.0000e-06
    0 : |proj(x-g)-x|    = 6.2729e+01 <= tolG          = 1.0000e-01
    0 : |proj(x-g)-x|    = 6.2729e+01 <= 1e3*eps       = 1.0000e-02
    0 : maxIter   =      50    <= iter          =      5
    ------------------------- DONE! -------------------------
    /home/vsts/work/1/s/examples/10-pgi/plot_inv_1_PGI_Linear_1D_joint_WithRelationships.py:301: UserWarning:

    marker is redundantly defined by the 'marker' keyword argument and the fmt string "b.-" (-> marker='.'). The keyword argument will take precedence.

    /home/vsts/work/1/s/examples/10-pgi/plot_inv_1_PGI_Linear_1D_joint_WithRelationships.py:308: UserWarning:

    marker is redundantly defined by the 'marker' keyword argument and the fmt string "r.-" (-> marker='.'). The keyword argument will take precedence.

    /home/vsts/work/1/s/examples/10-pgi/plot_inv_1_PGI_Linear_1D_joint_WithRelationships.py:346: UserWarning:

    marker is redundantly defined by the 'marker' keyword argument and the fmt string "b.-" (-> marker='.'). The keyword argument will take precedence.

    /home/vsts/work/1/s/examples/10-pgi/plot_inv_1_PGI_Linear_1D_joint_WithRelationships.py:353: UserWarning:

    marker is redundantly defined by the 'marker' keyword argument and the fmt string "r.-" (-> marker='.'). The keyword argument will take precedence.

    /home/vsts/work/1/s/examples/10-pgi/plot_inv_1_PGI_Linear_1D_joint_WithRelationships.py:360: UserWarning:

    The following kwargs were not used by contour: 'label'

    /home/vsts/work/1/s/examples/10-pgi/plot_inv_1_PGI_Linear_1D_joint_WithRelationships.py:368: UserWarning:

    The following kwargs were not used by contour: 'label'

    /home/vsts/work/1/s/examples/10-pgi/plot_inv_1_PGI_Linear_1D_joint_WithRelationships.py:402: UserWarning:

    marker is redundantly defined by the 'marker' keyword argument and the fmt string "b.-" (-> marker='.'). The keyword argument will take precedence.

    /home/vsts/work/1/s/examples/10-pgi/plot_inv_1_PGI_Linear_1D_joint_WithRelationships.py:409: UserWarning:

    marker is redundantly defined by the 'marker' keyword argument and the fmt string "r.-" (-> marker='.'). The keyword argument will take precedence.







|

.. code-block:: default


    import discretize as Mesh
    import matplotlib.pyplot as plt
    import numpy as np
    from SimPEG import (
        data_misfit,
        directives,
        inverse_problem,
        inversion,
        maps,
        optimization,
        regularization,
        simulation,
        utils,
    )

    # Random seed for reproductibility
    np.random.seed(1)
    # Mesh
    N = 100
    mesh = Mesh.TensorMesh([N])

    # Survey design parameters
    nk = 30
    jk = np.linspace(1.0, 59.0, nk)
    p = -0.25
    q = 0.25


    # Physics
    def g(k):
        return np.exp(p * jk[k] * mesh.cell_centers_x) * np.cos(
            np.pi * q * jk[k] * mesh.cell_centers_x
        )


    G = np.empty((nk, mesh.nC))

    for i in range(nk):
        G[i, :] = g(i)

    m0 = np.zeros(mesh.nC)
    m0[20:41] = np.linspace(0.0, 1.0, 21)
    m0[41:57] = np.linspace(-1, 0.0, 16)

    poly0 = maps.PolynomialPetroClusterMap(coeffyx=np.r_[0.0, -4.0, 4.0])
    poly1 = maps.PolynomialPetroClusterMap(coeffyx=np.r_[-0.0, 3.0, 6.0, 6.0])
    poly0_inverse = maps.PolynomialPetroClusterMap(coeffyx=-np.r_[0.0, -4.0, 4.0])
    poly1_inverse = maps.PolynomialPetroClusterMap(coeffyx=-np.r_[0.0, 3.0, 6.0, 6.0])
    cluster_mapping = [maps.IdentityMap(), poly0_inverse, poly1_inverse]

    m1 = np.zeros(100)
    m1[20:41] = 1.0 + (poly0 * np.vstack([m0[20:41], m1[20:41]]).T)[:, 1]
    m1[41:57] = -1.0 + (poly1 * np.vstack([m0[41:57], m1[41:57]]).T)[:, 1]

    model2d = np.vstack([m0, m1]).T
    m = utils.mkvc(model2d)

    clfmapping = utils.GaussianMixtureWithNonlinearRelationships(
        mesh=mesh,
        n_components=3,
        covariance_type="full",
        tol=1e-8,
        reg_covar=1e-3,
        max_iter=1000,
        n_init=100,
        init_params="kmeans",
        random_state=None,
        warm_start=False,
        means_init=np.array(
            [
                [0, 0],
                [m0[20:41].mean(), m1[20:41].mean()],
                [m0[41:57].mean(), m1[41:57].mean()],
            ]
        ),
        verbose=0,
        verbose_interval=10,
        cluster_mapping=cluster_mapping,
    )
    clfmapping = clfmapping.fit(model2d)

    clfnomapping = utils.WeightedGaussianMixture(
        mesh=mesh,
        n_components=3,
        covariance_type="full",
        tol=1e-8,
        reg_covar=1e-3,
        max_iter=1000,
        n_init=100,
        init_params="kmeans",
        random_state=None,
        warm_start=False,
        verbose=0,
        verbose_interval=10,
    )
    clfnomapping = clfnomapping.fit(model2d)

    wires = maps.Wires(("m1", mesh.nC), ("m2", mesh.nC))

    relatrive_error = 0.01
    noise_floor = 0.0

    prob1 = simulation.LinearSimulation(mesh, G=G, model_map=wires.m1)
    survey1 = prob1.make_synthetic_data(
        m, relative_error=relatrive_error, noise_floor=noise_floor, add_noise=True
    )

    prob2 = simulation.LinearSimulation(mesh, G=G, model_map=wires.m2)
    survey2 = prob2.make_synthetic_data(
        m, relative_error=relatrive_error, noise_floor=noise_floor, add_noise=True
    )


    dmis1 = data_misfit.L2DataMisfit(simulation=prob1, data=survey1)
    dmis2 = data_misfit.L2DataMisfit(simulation=prob2, data=survey2)
    dmis = dmis1 + dmis2
    minit = np.zeros_like(m)

    # Distance weighting
    wr1 = np.sum(prob1.G**2.0, axis=0) ** 0.5 / mesh.cell_volumes
    wr1 = wr1 / np.max(wr1)
    wr2 = np.sum(prob2.G**2.0, axis=0) ** 0.5 / mesh.cell_volumes
    wr2 = wr2 / np.max(wr2)

    reg_simple = regularization.PGI(
        mesh=mesh,
        gmmref=clfmapping,
        gmm=clfmapping,
        approx_gradient=True,
        wiresmap=wires,
        non_linear_relationships=True,
        weights_list=[wr1, wr2],
    )

    opt = optimization.ProjectedGNCG(
        maxIter=50,
        tolX=1e-6,
        maxIterCG=100,
        tolCG=1e-3,
        lower=-10,
        upper=10,
    )

    invProb = inverse_problem.BaseInvProblem(dmis, reg_simple, opt)

    # directives
    scales = directives.ScalingMultipleDataMisfits_ByEig(
        chi0_ratio=np.r_[1.0, 1.0], verbose=True, n_pw_iter=10
    )
    scaling_schedule = directives.JointScalingSchedule(verbose=True)
    alpha0_ratio = np.r_[1e6, 1e4, 1, 1]
    alphas = directives.AlphasSmoothEstimate_ByEig(
        alpha0_ratio=alpha0_ratio, n_pw_iter=10, verbose=True
    )
    beta = directives.BetaEstimate_ByEig(beta0_ratio=1e-5, n_pw_iter=10)
    betaIt = directives.PGI_BetaAlphaSchedule(
        verbose=True,
        coolingFactor=2.0,
        progress=0.2,
    )
    targets = directives.MultiTargetMisfits(verbose=True)
    petrodir = directives.PGI_UpdateParameters(update_gmm=False)

    # Setup Inversion
    inv = inversion.BaseInversion(
        invProb,
        directiveList=[alphas, scales, beta, petrodir, targets, betaIt, scaling_schedule],
    )

    mcluster_map = inv.run(minit)

    # Inversion with no nonlinear mapping
    reg_simple_no_map = regularization.PGI(
        mesh=mesh,
        gmmref=clfnomapping,
        gmm=clfnomapping,
        approx_gradient=True,
        wiresmap=wires,
        non_linear_relationships=False,
        weights_list=[wr1, wr2],
    )

    opt = optimization.ProjectedGNCG(
        maxIter=50,
        tolX=1e-6,
        maxIterCG=100,
        tolCG=1e-3,
        lower=-10,
        upper=10,
    )

    invProb = inverse_problem.BaseInvProblem(dmis, reg_simple_no_map, opt)

    # directives
    scales = directives.ScalingMultipleDataMisfits_ByEig(
        chi0_ratio=np.r_[1.0, 1.0], verbose=True, n_pw_iter=10
    )
    scaling_schedule = directives.JointScalingSchedule(verbose=True)
    alpha0_ratio = np.r_[100.0 * np.ones(2), 1, 1]
    alphas = directives.AlphasSmoothEstimate_ByEig(
        alpha0_ratio=alpha0_ratio, n_pw_iter=10, verbose=True
    )
    beta = directives.BetaEstimate_ByEig(beta0_ratio=1e-5, n_pw_iter=10)
    betaIt = directives.PGI_BetaAlphaSchedule(
        verbose=True,
        coolingFactor=2.0,
        progress=0.2,
    )
    targets = directives.MultiTargetMisfits(
        chiSmall=1.0, TriggerSmall=True, TriggerTheta=False, verbose=True
    )
    petrodir = directives.PGI_UpdateParameters(update_gmm=False)

    # Setup Inversion
    inv = inversion.BaseInversion(
        invProb,
        directiveList=[alphas, scales, beta, petrodir, targets, betaIt, scaling_schedule],
    )

    mcluster_no_map = inv.run(minit)

    # WeightedLeastSquares Inversion

    reg1 = regularization.WeightedLeastSquares(
        mesh, alpha_s=1.0, alpha_x=1.0, mapping=wires.m1
    )
    reg1.cell_weights = wr1
    reg2 = regularization.WeightedLeastSquares(
        mesh, alpha_s=1.0, alpha_x=1.0, mapping=wires.m2
    )
    reg2.cell_weights = wr2
    reg = reg1 + reg2

    opt = optimization.ProjectedGNCG(
        maxIter=50,
        tolX=1e-6,
        maxIterCG=100,
        tolCG=1e-3,
        lower=-10,
        upper=10,
    )

    invProb = inverse_problem.BaseInvProblem(dmis, reg, opt)

    # directives
    alpha0_ratio = np.r_[1, 1, 1, 1]
    alphas = directives.AlphasSmoothEstimate_ByEig(
        alpha0_ratio=alpha0_ratio, n_pw_iter=10, verbose=True
    )
    scales = directives.ScalingMultipleDataMisfits_ByEig(
        chi0_ratio=np.r_[1.0, 1.0], verbose=True, n_pw_iter=10
    )
    scaling_schedule = directives.JointScalingSchedule(verbose=True)
    beta = directives.BetaEstimate_ByEig(beta0_ratio=1e-5, n_pw_iter=10)
    beta_schedule = directives.BetaSchedule(coolingFactor=5.0, coolingRate=1)
    targets = directives.MultiTargetMisfits(
        TriggerSmall=False,
        verbose=True,
    )

    # Setup Inversion
    inv = inversion.BaseInversion(
        invProb,
        directiveList=[alphas, scales, beta, targets, beta_schedule, scaling_schedule],
    )

    mtik = inv.run(minit)


    # Final Plot
    fig, axes = plt.subplots(3, 4, figsize=(25, 15))
    axes = axes.reshape(12)
    left, width = 0.25, 0.5
    bottom, height = 0.25, 0.5
    right = left + width
    top = bottom + height

    axes[0].set_axis_off()
    axes[0].text(
        0.5 * (left + right),
        0.5 * (bottom + top),
        ("Using true nonlinear\npetrophysical relationships"),
        horizontalalignment="center",
        verticalalignment="center",
        fontsize=20,
        color="black",
        transform=axes[0].transAxes,
    )

    axes[1].plot(mesh.cell_centers_x, wires.m1 * mcluster_map, "b.-", ms=5, marker="v")
    axes[1].plot(mesh.cell_centers_x, wires.m1 * m, "k--")
    axes[1].set_title("Problem 1")
    axes[1].legend(["Recovered Model", "True Model"], loc=1)
    axes[1].set_xlabel("X")
    axes[1].set_ylabel("Property 1")

    axes[2].plot(mesh.cell_centers_x, wires.m2 * mcluster_map, "r.-", ms=5, marker="v")
    axes[2].plot(mesh.cell_centers_x, wires.m2 * m, "k--")
    axes[2].set_title("Problem 2")
    axes[2].legend(["Recovered Model", "True Model"], loc=1)
    axes[2].set_xlabel("X")
    axes[2].set_ylabel("Property 2")

    x, y = np.mgrid[-1:1:0.01, -4:2:0.01]
    pos = np.empty(x.shape + (2,))
    pos[:, :, 0] = x
    pos[:, :, 1] = y
    CS = axes[3].contour(
        x,
        y,
        np.exp(clfmapping.score_samples(pos.reshape(-1, 2)).reshape(x.shape)),
        100,
        alpha=0.25,
        cmap="viridis",
    )
    axes[3].scatter(wires.m1 * mcluster_map, wires.m2 * mcluster_map, marker="v")
    axes[3].set_title("Petrophysical Distribution")
    CS.collections[0].set_label("")
    axes[3].legend(["True Petrophysical Distribution", "Recovered model crossplot"])
    axes[3].set_xlabel("Property 1")
    axes[3].set_ylabel("Property 2")

    axes[4].set_axis_off()
    axes[4].text(
        0.5 * (left + right),
        0.5 * (bottom + top),
        ("Using a pure\nGaussian distribution"),
        horizontalalignment="center",
        verticalalignment="center",
        fontsize=20,
        color="black",
        transform=axes[4].transAxes,
    )

    axes[5].plot(mesh.cell_centers_x, wires.m1 * mcluster_no_map, "b.-", ms=5, marker="v")
    axes[5].plot(mesh.cell_centers_x, wires.m1 * m, "k--")
    axes[5].set_title("Problem 1")
    axes[5].legend(["Recovered Model", "True Model"], loc=1)
    axes[5].set_xlabel("X")
    axes[5].set_ylabel("Property 1")

    axes[6].plot(mesh.cell_centers_x, wires.m2 * mcluster_no_map, "r.-", ms=5, marker="v")
    axes[6].plot(mesh.cell_centers_x, wires.m2 * m, "k--")
    axes[6].set_title("Problem 2")
    axes[6].legend(["Recovered Model", "True Model"], loc=1)
    axes[6].set_xlabel("X")
    axes[6].set_ylabel("Property 2")

    CSF = axes[7].contour(
        x,
        y,
        np.exp(clfmapping.score_samples(pos.reshape(-1, 2)).reshape(x.shape)),
        100,
        alpha=0.5,
        label="True Petro. Distribution",
    )
    CS = axes[7].contour(
        x,
        y,
        np.exp(clfnomapping.score_samples(pos.reshape(-1, 2)).reshape(x.shape)),
        500,
        cmap="viridis",
        linestyles="--",
        label="Modeled Petro. Distribution",
    )
    axes[7].scatter(
        wires.m1 * mcluster_no_map,
        wires.m2 * mcluster_no_map,
        marker="v",
        label="Recovered model crossplot",
    )
    axes[7].set_title("Petrophysical Distribution")
    axes[7].legend()
    axes[7].set_xlabel("Property 1")
    axes[7].set_ylabel("Property 2")

    # Tikonov

    axes[8].set_axis_off()
    axes[8].text(
        0.5 * (left + right),
        0.5 * (bottom + top),
        ("Least-Squares\n~Using a single cluster"),
        horizontalalignment="center",
        verticalalignment="center",
        fontsize=20,
        color="black",
        transform=axes[8].transAxes,
    )

    axes[9].plot(mesh.cell_centers_x, wires.m1 * mtik, "b.-", ms=5, marker="v")
    axes[9].plot(mesh.cell_centers_x, wires.m1 * m, "k--")
    axes[9].set_title("Problem 1")
    axes[9].legend(["Recovered Model", "True Model"], loc=1)
    axes[9].set_xlabel("X")
    axes[9].set_ylabel("Property 1")

    axes[10].plot(mesh.cell_centers_x, wires.m2 * mtik, "r.-", ms=5, marker="v")
    axes[10].plot(mesh.cell_centers_x, wires.m2 * m, "k--")
    axes[10].set_title("Problem 2")
    axes[10].legend(["Recovered Model", "True Model"], loc=1)
    axes[10].set_xlabel("X")
    axes[10].set_ylabel("Property 2")

    CS = axes[11].contour(
        x,
        y,
        np.exp(clfmapping.score_samples(pos.reshape(-1, 2)).reshape(x.shape)),
        100,
        alpha=0.25,
        cmap="viridis",
    )
    axes[11].scatter(wires.m1 * mtik, wires.m2 * mtik, marker="v")
    axes[11].set_title("Petro Distribution")
    CS.collections[0].set_label("")
    axes[11].legend(["True Petro Distribution", "Recovered model crossplot"])
    axes[11].set_xlabel("Property 1")
    axes[11].set_ylabel("Property 2")
    plt.subplots_adjust(wspace=0.3, hspace=0.3, top=0.85)
    plt.show()


.. rst-class:: sphx-glr-timing

   **Total running time of the script:** (0 minutes 47.091 seconds)

**Estimated memory usage:**  8 MB


.. _sphx_glr_download_content_examples_10-pgi_plot_inv_1_PGI_Linear_1D_joint_WithRelationships.py:

.. only:: html

  .. container:: sphx-glr-footer sphx-glr-footer-example




    .. container:: sphx-glr-download sphx-glr-download-python

      :download:`Download Python source code: plot_inv_1_PGI_Linear_1D_joint_WithRelationships.py <plot_inv_1_PGI_Linear_1D_joint_WithRelationships.py>`

    .. container:: sphx-glr-download sphx-glr-download-jupyter

      :download:`Download Jupyter notebook: plot_inv_1_PGI_Linear_1D_joint_WithRelationships.ipynb <plot_inv_1_PGI_Linear_1D_joint_WithRelationships.ipynb>`


.. only:: html

 .. rst-class:: sphx-glr-signature

    `Gallery generated by Sphinx-Gallery <https://sphinx-gallery.github.io>`_