.. 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 ` 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 Running inversion with SimPEG v0.23.0 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: [np.float64(3.5076401722247157), np.float64(0.0), np.float64(3.508490434878856e-06), np.float64(0.0)] Calculating the scaling parameter. Scale Multipliers: [0.09326744 0.90673256] Initial data misfit scales: [0.09326744 0.90673256] 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.87e+01 3.00e+05 0.00e+00 3.00e+05 1.41e+02 0 geophys. misfits: 1096.0 (target 30.0 [False]); 70.7 (target 30.0 [False]) | smallness misfit: 2986.0 (target: 200.0 [False]) Beta cooling evaluation: progress: [1096. 70.7]; minimum progress targets: [240000. 240000.] 1 1.87e+01 1.66e+02 4.16e+01 9.44e+02 8.82e+01 0 geophys. misfits: 510.0 (target 30.0 [False]); 18.9 (target 30.0 [True]) | smallness misfit: 1373.3 (target: 200.0 [False]) Beta cooling evaluation: progress: [510. 18.9]; minimum progress targets: [876.8 56.5] Updating scaling for data misfits by 1.5859329981033554 New scales: [0.14025139 0.85974861] 2 1.87e+01 8.78e+01 4.05e+01 8.44e+02 8.95e+01 0 Skip BFGS geophys. misfits: 291.8 (target 30.0 [False]); 18.6 (target 30.0 [True]) | smallness misfit: 1269.4 (target: 200.0 [False]) Beta cooling evaluation: progress: [291.8 18.6]; minimum progress targets: [408. 30.] Updating scaling for data misfits by 1.6128264996374915 New scales: [0.208298 0.791702] 3 1.87e+01 7.55e+01 4.18e+01 8.56e+02 7.51e+01 0 Skip BFGS geophys. misfits: 172.5 (target 30.0 [False]); 18.7 (target 30.0 [True]) | smallness misfit: 1200.0 (target: 200.0 [False]) Beta cooling evaluation: progress: [172.5 18.7]; minimum progress targets: [233.5 30. ] Updating scaling for data misfits by 1.6048541987288436 New scales: [0.29688358 0.70311642] 4 1.87e+01 6.43e+01 4.28e+01 8.65e+02 7.47e+01 0 geophys. misfits: 111.6 (target 30.0 [False]); 19.2 (target 30.0 [True]) | smallness misfit: 1148.5 (target: 200.0 [False]) Beta cooling evaluation: progress: [111.6 19.2]; minimum progress targets: [138. 30.] Updating scaling for data misfits by 1.5624977847730535 New scales: [0.39749904 0.60250096] 5 1.87e+01 5.59e+01 4.36e+01 8.71e+02 8.06e+01 0 Skip BFGS geophys. misfits: 80.1 (target 30.0 [False]); 20.2 (target 30.0 [True]) | smallness misfit: 1106.0 (target: 200.0 [False]) Beta cooling evaluation: progress: [80.1 20.2]; minimum progress targets: [89.3 30. ] Updating scaling for data misfits by 1.4865180376781295 New scales: [0.49513509 0.50486491] 6 1.87e+01 4.99e+01 4.41e+01 8.75e+02 8.43e+01 0 Skip BFGS geophys. misfits: 63.7 (target 30.0 [False]); 21.7 (target 30.0 [True]) | smallness misfit: 1066.6 (target: 200.0 [False]) Beta cooling evaluation: progress: [63.7 21.7]; minimum progress targets: [64.1 30. ] Updating scaling for data misfits by 1.3800352190916845 New scales: [0.5750899 0.4249101] 7 1.87e+01 4.58e+01 4.45e+01 8.78e+02 8.33e+01 0 Skip BFGS geophys. misfits: 55.1 (target 30.0 [False]); 24.1 (target 30.0 [True]) | smallness misfit: 1029.5 (target: 200.0 [False]) Beta cooling evaluation: progress: [55.1 24.1]; minimum progress targets: [50.9 30. ] Decreasing beta to counter data misfit decrase plateau. Updating scaling for data misfits by 1.2469048359310606 New scales: [0.62792216 0.37207784] 8 9.35e+00 4.36e+01 4.47e+01 4.61e+02 9.42e+01 0 Skip BFGS geophys. misfits: 30.1 (target 30.0 [False]); 17.5 (target 30.0 [True]) | smallness misfit: 1105.5 (target: 200.0 [False]) Beta cooling evaluation: progress: [30.1 17.5]; minimum progress targets: [44.1 30. ] Updating scaling for data misfits by 1.7174836259931823 New scales: [0.74348727 0.25651273] 9 9.35e+00 2.68e+01 4.59e+01 4.56e+02 7.63e+01 0 geophys. misfits: 27.5 (target 30.0 [True]); 20.4 (target 30.0 [True]) | smallness misfit: 1044.1 (target: 200.0 [False]) Beta cooling evaluation: progress: [27.5 20.4]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 1.2786296516513311 10 9.35e+00 2.57e+01 4.67e+01 4.62e+02 7.63e+01 0 Skip BFGS geophys. misfits: 27.4 (target 30.0 [True]); 23.0 (target 30.0 [True]) | smallness misfit: 973.7 (target: 200.0 [False]) Beta cooling evaluation: progress: [27.4 23. ]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 1.5327952864032601 11 9.35e+00 2.63e+01 4.71e+01 4.67e+02 7.60e+01 0 geophys. misfits: 27.3 (target 30.0 [True]); 24.9 (target 30.0 [True]) | smallness misfit: 927.5 (target: 200.0 [False]) Beta cooling evaluation: progress: [27.3 24.9]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 1.7667402301739434 12 9.35e+00 2.67e+01 4.75e+01 4.71e+02 7.71e+01 0 Skip BFGS geophys. misfits: 27.2 (target 30.0 [True]); 27.2 (target 30.0 [True]) | smallness misfit: 884.1 (target: 200.0 [False]) Beta cooling evaluation: progress: [27.2 27.2]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 1.9511303520110128 13 9.35e+00 2.72e+01 4.78e+01 4.75e+02 7.54e+01 0 geophys. misfits: 27.4 (target 30.0 [True]); 29.5 (target 30.0 [True]) | smallness misfit: 856.6 (target: 200.0 [False]) Beta cooling evaluation: progress: [27.4 29.5]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 2.063663142327159 14 9.35e+00 2.79e+01 4.80e+01 4.76e+02 7.55e+01 0 Skip BFGS geophys. misfits: 27.1 (target 30.0 [True]); 30.6 (target 30.0 [False]) | smallness misfit: 841.2 (target: 200.0 [False]) Beta cooling evaluation: progress: [27.1 30.6]; minimum progress targets: [30. 30.] Decreasing beta to counter data misfit increase. Updating scaling for data misfits by 1.1050286910055975 New scales: [0.72398232 0.27601768] 15 4.68e+00 2.81e+01 4.80e+01 2.52e+02 9.47e+01 0 geophys. misfits: 21.5 (target 30.0 [True]); 19.4 (target 30.0 [True]) | smallness misfit: 944.9 (target: 200.0 [False]) Beta cooling evaluation: progress: [21.5 19.4]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 3.0324546470037728 16 4.68e+00 2.09e+01 5.10e+01 2.59e+02 8.39e+01 0 geophys. misfits: 21.5 (target 30.0 [True]); 22.3 (target 30.0 [True]) | smallness misfit: 837.9 (target: 200.0 [False]) Beta cooling evaluation: progress: [21.5 22.3]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 4.149483792621317 17 4.68e+00 2.18e+01 5.27e+01 2.68e+02 5.66e+01 0 geophys. misfits: 21.5 (target 30.0 [True]); 23.9 (target 30.0 [True]) | smallness misfit: 766.0 (target: 200.0 [False]) Beta cooling evaluation: progress: [21.5 23.9]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 5.4994530769824985 18 4.68e+00 2.22e+01 5.49e+01 2.79e+02 8.52e+01 0 Skip BFGS geophys. misfits: 21.6 (target 30.0 [True]); 25.9 (target 30.0 [True]) | smallness misfit: 690.2 (target: 200.0 [False]) Beta cooling evaluation: progress: [21.6 25.9]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 6.997932378586193 19 4.68e+00 2.28e+01 5.69e+01 2.89e+02 9.58e+01 0 geophys. misfits: 21.8 (target 30.0 [True]); 29.8 (target 30.0 [True]) | smallness misfit: 648.9 (target: 200.0 [False]) Beta cooling evaluation: progress: [21.8 29.8]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 8.345249132782294 20 4.68e+00 2.40e+01 5.84e+01 2.97e+02 9.18e+01 0 geophys. misfits: 22.2 (target 30.0 [True]); 29.2 (target 30.0 [True]) | smallness misfit: 578.2 (target: 200.0 [False]) Beta cooling evaluation: progress: [22.2 29.2]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 9.909691752145934 21 4.68e+00 2.42e+01 6.01e+01 3.05e+02 9.66e+01 0 geophys. misfits: 22.9 (target 30.0 [True]); 35.7 (target 30.0 [False]) | smallness misfit: 487.2 (target: 200.0 [False]) Beta cooling evaluation: progress: [22.9 35.7]; minimum progress targets: [30. 30.] Decreasing beta to counter data misfit increase. Updating scaling for data misfits by 1.3078380663123905 New scales: [0.66728403 0.33271597] 22 2.34e+00 2.72e+01 5.92e+01 1.66e+02 1.07e+02 0 Skip BFGS geophys. misfits: 20.3 (target 30.0 [True]); 23.5 (target 30.0 [True]) | smallness misfit: 553.6 (target: 200.0 [False]) Beta cooling evaluation: progress: [20.3 23.5]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 13.643321046096773 23 2.34e+00 2.14e+01 6.54e+01 1.74e+02 7.84e+01 0 geophys. misfits: 20.6 (target 30.0 [True]); 21.0 (target 30.0 [True]) | smallness misfit: 509.5 (target: 200.0 [False]) Beta cooling evaluation: progress: [20.6 21. ]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 19.667094190799734 24 2.34e+00 2.07e+01 7.20e+01 1.89e+02 1.00e+02 0 geophys. misfits: 20.6 (target 30.0 [True]); 19.2 (target 30.0 [True]) | smallness misfit: 468.0 (target: 200.0 [False]) Beta cooling evaluation: progress: [20.6 19.2]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 29.74263976375323 25 2.34e+00 2.01e+01 8.22e+01 2.12e+02 1.14e+02 0 geophys. misfits: 20.6 (target 30.0 [True]); 18.2 (target 30.0 [True]) | smallness misfit: 426.3 (target: 200.0 [False]) Beta cooling evaluation: progress: [20.6 18.2]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 46.136205946409 26 2.34e+00 1.98e+01 9.69e+01 2.46e+02 1.24e+02 0 geophys. misfits: 21.6 (target 30.0 [True]); 18.0 (target 30.0 [True]) | smallness misfit: 384.2 (target: 200.0 [False]) Beta cooling evaluation: progress: [21.6 18. ]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 70.45051841767439 27 2.34e+00 2.04e+01 1.16e+02 2.92e+02 1.29e+02 0 geophys. misfits: 22.5 (target 30.0 [True]); 19.0 (target 30.0 [True]) | smallness misfit: 342.8 (target: 200.0 [False]) Beta cooling evaluation: progress: [22.5 19. ]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 102.63896119771891 28 2.34e+00 2.13e+01 1.38e+02 3.44e+02 1.31e+02 0 geophys. misfits: 22.4 (target 30.0 [True]); 25.0 (target 30.0 [True]) | smallness misfit: 332.4 (target: 200.0 [False]) Beta cooling evaluation: progress: [22.4 25. ]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 130.30166640202273 29 2.34e+00 2.33e+01 1.57e+02 3.89e+02 1.31e+02 0 geophys. misfits: 27.5 (target 30.0 [True]); 28.6 (target 30.0 [True]) | smallness misfit: 272.7 (target: 200.0 [False]) Beta cooling evaluation: progress: [27.5 28.6]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 139.33865961698652 30 2.34e+00 2.79e+01 1.57e+02 3.95e+02 1.34e+02 0 geophys. misfits: 25.4 (target 30.0 [True]); 27.1 (target 30.0 [True]) | smallness misfit: 260.1 (target: 200.0 [False]) Beta cooling evaluation: progress: [25.4 27.1]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 159.4082639968159 31 2.34e+00 2.60e+01 1.61e+02 4.02e+02 1.31e+02 1 geophys. misfits: 24.3 (target 30.0 [True]); 25.9 (target 30.0 [True]) | smallness misfit: 231.7 (target: 200.0 [False]) Beta cooling evaluation: progress: [24.3 25.9]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 190.80460809699628 32 2.34e+00 2.48e+01 1.68e+02 4.17e+02 1.32e+02 1 geophys. misfits: 28.7 (target 30.0 [True]); 50.5 (target 30.0 [False]) | smallness misfit: 174.4 (target: 200.0 [True]) Beta cooling evaluation: progress: [28.7 50.5]; minimum progress targets: [30. 30.] Decreasing beta to counter data misfit increase. Updating scaling for data misfits by 1.0459707133398461 New scales: [0.65723155 0.34276845] 33 1.17e+00 3.62e+01 1.55e+02 2.17e+02 1.25e+02 0 geophys. misfits: 23.8 (target 30.0 [True]); 59.5 (target 30.0 [False]) | smallness misfit: 153.9 (target: 200.0 [True]) Beta cooling evaluation: progress: [23.8 59.5]; minimum progress targets: [30. 40.4] Decreasing beta to counter data misfit increase. Updating scaling for data misfits by 1.2596773070109566 New scales: [0.60351331 0.39648669] 34 5.84e-01 3.80e+01 1.43e+02 1.22e+02 1.10e+02 0 geophys. misfits: 22.2 (target 30.0 [True]); 25.6 (target 30.0 [True]) | smallness misfit: 178.8 (target: 200.0 [True]) All targets have been reached Beta cooling evaluation: progress: [22.2 25.6]; minimum progress targets: [30. 47.6] Warming alpha_pgi to favor clustering: 240.72513045838488 ------------------------- STOP! ------------------------- 1 : |fc-fOld| = 0.0000e+00 <= tolF*(1+|f0|) = 3.0000e+04 0 : |xc-x_last| = 7.9996e-01 <= tolX*(1+|x0|) = 1.0000e-06 0 : |proj(x-g)-x| = 1.0991e+02 <= tolG = 1.0000e-01 0 : |proj(x-g)-x| = 1.0991e+02 <= 1e3*eps = 1.0000e-02 0 : maxIter = 50 <= iter = 35 ------------------------- DONE! ------------------------- Running inversion with SimPEG v0.23.0 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: [np.float64(0.00039124967802256705), np.float64(0.0), np.float64(3.4991079985601646e-06), np.float64(0.0)] Calculating the scaling parameter. Scale Multipliers: [0.09326744 0.90673256] Initial data misfit scales: [0.09326744 0.90673256] 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.93e+03 3.00e+05 0.00e+00 3.00e+05 1.41e+02 0 geophys. misfits: 89877.2 (target 30.0 [False]); 62289.4 (target 30.0 [False]) | smallness misfit: 294.2 (target: 200.0 [False]) Beta cooling evaluation: progress: [89877.2 62289.4]; minimum progress targets: [240000. 240000.] 1 1.93e+03 6.49e+04 6.93e-01 6.62e+04 9.06e+01 0 geophys. misfits: 675.2 (target 30.0 [False]); 20.1 (target 30.0 [True]) | smallness misfit: 96.5 (target: 200.0 [True]) Beta cooling evaluation: progress: [675.2 20.1]; minimum progress targets: [71901.8 49831.5] Updating scaling for data misfits by 1.493484799701024 New scales: [0.13316447 0.86683553] 2 1.93e+03 1.07e+02 2.53e-01 5.95e+02 9.57e+01 0 Skip BFGS geophys. misfits: 94.0 (target 30.0 [False]); 18.9 (target 30.0 [True]) | smallness misfit: 51.3 (target: 200.0 [True]) Beta cooling evaluation: progress: [94. 18.9]; minimum progress targets: [540.1 30. ] Updating scaling for data misfits by 1.5900801673882476 New scales: [0.19631611 0.80368389] 3 1.93e+03 3.36e+01 1.46e-01 3.16e+02 1.06e+02 0 Skip BFGS geophys. misfits: 59.7 (target 30.0 [False]); 18.9 (target 30.0 [True]) | smallness misfit: 47.9 (target: 200.0 [True]) Beta cooling evaluation: progress: [59.7 18.9]; minimum progress targets: [75.2 30. ] Updating scaling for data misfits by 1.588985131008552 New scales: [0.27961255 0.72038745] 4 1.93e+03 3.03e+01 1.31e-01 2.84e+02 7.87e+01 0 geophys. misfits: 44.2 (target 30.0 [False]); 19.2 (target 30.0 [True]) | smallness misfit: 48.7 (target: 200.0 [True]) Beta cooling evaluation: progress: [44.2 19.2]; minimum progress targets: [47.8 30. ] Updating scaling for data misfits by 1.5626108355934907 New scales: [0.3775345 0.6224655] 5 1.93e+03 2.86e+01 1.33e-01 2.85e+02 8.45e+01 0 Skip BFGS geophys. misfits: 36.3 (target 30.0 [False]); 20.4 (target 30.0 [True]) | smallness misfit: 49.1 (target: 200.0 [True]) Beta cooling evaluation: progress: [36.3 20.4]; minimum progress targets: [35.4 30. ] Decreasing beta to counter data misfit decrase plateau. Updating scaling for data misfits by 1.4725735801347535 New scales: [0.47177637 0.52822363] 6 9.64e+02 2.79e+01 1.34e-01 1.57e+02 1.05e+02 0 geophys. misfits: 25.9 (target 30.0 [True]); 17.3 (target 30.0 [True]) | smallness misfit: 51.1 (target: 200.0 [True]) All targets have been reached Beta cooling evaluation: progress: [25.9 17.3]; minimum progress targets: [30. 30.] Warming alpha_pgi to favor clustering: 1.4442529544627405 ------------------------- STOP! ------------------------- 1 : |fc-fOld| = 0.0000e+00 <= tolF*(1+|f0|) = 3.0000e+04 0 : |xc-x_last| = 1.5816e-01 <= tolX*(1+|x0|) = 1.0000e-06 0 : |proj(x-g)-x| = 1.0514e+02 <= tolG = 1.0000e-01 0 : |proj(x-g)-x| = 1.0514e+02 <= 1e3*eps = 1.0000e-02 0 : maxIter = 50 <= iter = 7 ------------------------- DONE! ------------------------- Running inversion with SimPEG v0.23.0 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: [np.float64(3.510762811221696e-05), np.float64(0.0), np.float64(4.420996595587482e-05), np.float64(0.0)] Calculating the scaling parameter. Scale Multipliers: [0.09326744 0.90673256] /home/vsts/work/1/s/simpeg/directives/directives.py:339: UserWarning: There is no PGI regularization. Smallness target is turned off (TriggerSmall flag) Initial data misfit scales: [0.09326744 0.90673256] 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 3.00e+05 0.00e+00 3.00e+05 1.41e+02 0 geophys. misfits: 56827.8 (target 30.0 [False]); 34909.5 (target 30.0 [False]) 1 2.03e+05 3.70e+04 4.35e-02 4.58e+04 1.38e+02 0 geophys. misfits: 8507.7 (target 30.0 [False]); 3707.9 (target 30.0 [False]) 2 4.07e+04 4.16e+03 1.07e-01 8.53e+03 1.30e+02 0 Skip BFGS geophys. misfits: 582.4 (target 30.0 [False]); 239.5 (target 30.0 [False]) 3 8.14e+03 2.71e+02 1.41e-01 1.42e+03 1.04e+02 0 Skip BFGS geophys. misfits: 50.8 (target 30.0 [False]); 27.3 (target 30.0 [True]) Updating scaling for data misfits by 1.0986468594984395 New scales: [0.10153382 0.89846618] 4 1.63e+03 2.97e+01 1.52e-01 2.76e+02 8.03e+01 0 Skip BFGS geophys. misfits: 20.7 (target 30.0 [True]); 12.6 (target 30.0 [True]) All targets have been reached ------------------------- STOP! ------------------------- 1 : |fc-fOld| = 0.0000e+00 <= tolF*(1+|f0|) = 3.0000e+04 0 : |xc-x_last| = 4.3399e-01 <= tolX*(1+|x0|) = 1.0000e-06 0 : |proj(x-g)-x| = 8.0315e+01 <= tolG = 1.0000e-01 0 : |proj(x-g)-x| = 8.0315e+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:329: MatplotlibDeprecationWarning: The collections attribute was deprecated in Matplotlib 3.8 and will be removed in 3.10. /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. /home/vsts/work/1/s/examples/10-pgi/plot_inv_1_PGI_Linear_1D_joint_WithRelationships.py:426: MatplotlibDeprecationWarning: The collections attribute was deprecated in Matplotlib 3.8 and will be removed in 3.10. | .. code-block:: Python 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, weights={"cell_weights": wr1} ) reg2 = regularization.WeightedLeastSquares( mesh, alpha_s=1.0, alpha_x=1.0, mapping=wires.m2, weights={"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 56.858 seconds) **Estimated memory usage:** 288 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-jupyter :download:`Download Jupyter notebook: plot_inv_1_PGI_Linear_1D_joint_WithRelationships.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_inv_1_PGI_Linear_1D_joint_WithRelationships.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_inv_1_PGI_Linear_1D_joint_WithRelationships.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_