Faculty of Geo-Information Science and Earth Observation (ITC)
University of Twente, Langezijds Bld., Hallenweg 8, Enschede 7522 NH, the Netherlands.
TEL: + 61 (0)2 9850 8298
FAX: + 61 (0)2 9850 6904
e-mail: j.c.afonso@utwente.nl
SOFTWARE
Welcome!! Here you will find information about the latest versions of the "LitMod" codes and other software I have developed and/or contributed to.
If you are looking for the global model LithoRef18, please go to the bottom of the page. Find it also via the GPlates Portal
LitMod_2D and LitMod_T
LitMod_2D is an interactive software to perform multi-observable modelling of geophysical data for the whole lithospheric structure of the Earth and/or other terrestrial planets (in Windows and Linux). It was developed by J. C. Afonso but there are a few versions available from different groups. Currently, it is maintained mainly by the Group of Dynamics of the Lithosphere at the Institute of Earth Sciences J. Almera in Barcelona. You should get in touch with them if you are interested in using the latest versions of the code.
The main publication is:
- Kumar, A., Fernàndez, M., Jiménez-Munt, I., Torne, M., Vergés, J., & Afonso, J. C. (2020). LitMod2D_2.0: An improved integrated geophysical-petrological modeling tool for the physical interpretation of upper mantle anomalies. Geochemistry, Geophysics, Geosystems, 21, e2019GC008777. https://doi.org/10.1029/2019GC008777
Here are some relevant publications using LitMod_2D:
- Afonso, J.C., Fernandez, M., Ranalli, G., Griffin, W. L., Connolly, J. (2008). Combined geophysical-petrological modelling of the lithospheric-sublithospheric upper mantle: methodology and applications. Geochem. Geophys. Geosyst., 9, Q05008, doi:10.1029/2007GC001834. PDF
- Tunini, L., Jimenez-Munt, I., Fernandez, M., Verges, Villasenor, A. Melchiorre, M., Afonso, J.C., Geophysical-petrological model of the crust and upperm mantle in the India-Eurasia collision zone. Tectonics , 35, 1642-1669, doi:10.1002/2016TC004161 PDF
- Carballo, A., Fernandez, M., Jimenez-Munt, I., Torne, M., Verges, J., Melchiorre, M., Pedreira, D., Afonso, J.C., Garcia-Castellanos, D., Diaz, J., Villasenor, A., Pulgar, J.A., Quintana, L. (2015), From the Norht-Iberian Margin to the Alboran Basin: a lithosphere geo-transect crossing the Iberian Plate. Tectonophysics. , 663, 399-418. PDF
- Pedreira, D., Afonso, J.C. et al., (2015), Geophysical-petrological modeling of the lithosphere beneath the Cantabrian Mountains and North-Iberian margin: geodynamic implications. Lithos, 230, 46-68. PDF
- Fernandez, M., Afonso, J.C., Ranalli, G. (2010). The deep lithospheric structure of the Namibian volcanic margin Tectonophysics, 481, 68-81. PDF
- Alasonati-Tasarova, Z., Afonso, J.C., Bielik, M., Gotze, H-J. Hak, J. (2009). The lithospheric structure of the Western Carpathian-Pannonian Basin region based on the CELEBRATION 2000 seismic experiment and gravity modelling, Tectonophysics, 475, 454-469. PDF
Thermochemical model of western Tibet obtained with LitMod_2D by fitting surface heat flow, bouguer anomalies, geoid anomalies, elevation and seismic anomalies. f) anf g) show a comparison of predicted and observed P-wave anomalies. Results from Tunini et al. (2015).
LitMod_3D
LitMod_3D is the big brother of LitMod_2D and it's used in industry and academia. It is also fully interactive and it has been developed in collaboration with J. Fullea. The currently supported version run in Linux only. The latest version of the code and the user guide can be downloaded from the button below. LitMod_3D is distributed as a compressed tar file that contains all the required executables, sources, and shell script files. All sources are freely available. If you are interested in LitMod_3D, here are some relevant publications you may want to check:
- Fullea, J., Afonso, J.C., Connolly, J.A.D., Fernandez, M., Garcia-Castellanos, D., Zeyen, H. (2009). LitMod3D: an interactive 3D software to model the thermal, compositional, density, rheological, and seismological structure of the lithosphere and sublithospheric upper mantle. Geochem. Geophys. Geosyst, 10, Q08019, doi:10.1029/2009GC002391. PDF
- Fullea, J., Fernandez, M., Afonso, J.C., Verges, J., Zeyen, H. (2010). The structure and evolution of the lithosphere-asthenosphere boundary beneath the Atlanitc-Mediterranean Transition Region. Lithos, 120, 74-95. PDF
- Alasonati Tašárová, Z., J. Fullea, M. Bielik, and P. Środa (2016), Lithospheric structure of Central Europe: Puzzle pieces from Pannonian Basin to Trans-Europe an Suture Zone resolved by geophysical-petrological modeling, Tectonics, 35, doi:10.1002/2015TC003935. PDF
- Gradmann, S. Ebbing, J., Fullea, J. (2013) Integrated geophysical modelling of a lateral transition zone in the lithospheric mantle under Norway and Sweden, Geophys. J. Int., 194,1358–1373, https://doi.org/10.1093/gji/ggt213
A 3D rendering of the thermal structure of the Atlantic-Mediterranean region from the study of Fullea et al. (2010) using LitMod_3D.
LitMod_4INV
This is the "big daddy" of the LitMod suite. It is a fully nonlinear probabilistic inversion code(s) in spherical coordinates for the compositional and thermal structure of the lithosphere and upper mantle. It is capable of simultaneously inverting
gravity gradients, gravity anomalies, geoid height, surface heat flow, magnetotelluric data, receiver functions, surface-wave data, absolute elevation (including both static and dynamic effects) and body-wave data, together with
petrological information. It requires a decent cluster to run, although it is faster/cheaper when not all the above data sets are inverted together.
The code is fully functional, but under constant testing and development by J.C. Afonso's group at Macquarie University. Besides academics, LitMod_4INV is being used (and funded) by a number of industry partners and government agencies. Send me an email if you are interested in using it. Some relevant publications are:
- Jones, A.G., Afonso, J.C., Fullea, J., (2017), Geochemical and geophysical constrains on the dynamic topography of the Southern African Plateau., Geochem. Geophys. Geosyst., 18, 3556–3575, doi:10.1002/2017GC006908. PDF
- Qashqai, T.M., Afonso, J.C., Yang. Y. (2016), The crustal structure of the Arizona Transition Zone and southern Colorado Plateau from multi-observable probabilistic inversion..Geophys.Geochem.Geosys. , 17, doi:10.1002/2016GC006463. PDF
- Afonso, J.C., Rawlinson N., Yang, Y., Schutt, D., Jones, A.G., Fullea, J., Griffin, W.L. (2016), 3D multi-observable probabilistic inversion for the compositional and thermal structure of the lithosphere and upper mantle III: Thermochemical Tomography in the Western-Central US.J. Geophys. Res., 121, doi:10.1002/2016JB013049 PDF
Some results obtained with LitMod_4INV from Afonso et al. (2016) JGR.
- Guo, Z., Afonso, J.C., Qashqai, M., Yang, Y., Chen, J. (2016), Thermochemical structure of the North China Craton from multi-observable probabilistic inversion: extent and causes of cratonic lithosphere modification. Gondwana Res. ,37, 252-265. PDF
- Shan, B., Afonso, J.C., Yang, Y., Grose, C.J., Zheng, Y., Xiong, X. (2014), The compositional and thermal structure of the lithosphere and upper mantle beneath South China: Results from multi-observable probabilistic inversion. J. Geophys. Res., 119, 8417-8441, doi:10.1002/2014JB011412 PDF
- Afonso, J.C., Fullea J., Yang, Y., Connolly, J.A.D., Jones, A.G. (2013), 3D multi-observable probabilistic inversion for the compositional and thermal structure of the lithosphere and upper mantle II: General methodology and resolution analysis. J. Geophys. Res.,
doi:10.1002/jgrb.50123. PDF
- Afonso, J.C., Fullea J., Griffin, W.L., Yang, Y., Jones, A.G., Connolly, J.A.D., O'Reilly, S.Y. (2013), 3D multi-observable probabilistic inversion for the compositional and thermal structure of the lithosphere and upper mantle I: a priori information and geophysical observables. J. Geophys. Res., doi:10.1002/jgrb.50124. PDF
Efficient and general algorithm for implementing thermodynamic information into geodynamic and geophysical studies
Our approach is based on tensor rank decomposition methods, which transform the original multi-dimensional discrete information into a separated representation that contains significantly fewer terms, thus drastically reducing the amount of information to be stored in memory during a numerical simulation or geophysical inversion. Accordingly, the amount and resolution of the
thermodynamic information that can be used in a simulation or inversion increases substantially. Any thermodynamic software or database can be used to obtain the primary thermodynamic information.
Please refer to: Afonso, J.C., Zlotnik, S., Diez, P. (2015) “An efficient and general approach for implementing thermodynamic phase-equilibria information in geophysical and geodynamic studies”, Geochem, Geophys., Geosys., 16, doi:10.1002/2015GC006031
Density fields obtained from the energy minimization problem for water contents from 0 to 4% (from left to right) are shown in the first row. Same fields provided by the decomposition procedure with 100 terms are shown in the second row. Their absolute difference (full solution - separated solution) is shown in the third row.
ADR FEM code
This algorithm/code was presented in the paper:
- Oliveira, B., Afonso, J.C., Zlotnik, S. (2016), A Lagrangian-Eulerian finite element algorithm for advection-diffusion-reaction problems with phase change, Comput. Methods Appl. Mech. Engrg., 300, 375-401.
It is a matlab code that implements the algorithm described in the paper above, which is a very powerful and accurate algorithm to model a large number of problems that involve ADR with phase changes (e.g. metal casting, groundwater transport problems, magmatic systems, etc). The algorithm combines a Lagrangian formulation for the advection + reaction problem with the Eularian-based heat source method for the diffusion + phase change problem.
This algorithm is the basis for solving the energy equation of the multi-phase multi-component reactive transport problem tackled in Oliveira et al. (2017) and Oliveira et. al. (2018)... see below.
Contour plots of the three-body rotation problem. Advection (3 revolutions) and advection+diffusion (1 revolution) problems (by rows) for different interpolation schemes (by columns). The reference solutions are shown in the rightmost column.
Code for multi-phase multi-component reactive transport
This multi-algorithm method and code was presented in the paper:
- Oliveira, B., Afonso, J.C., Zlotnik, S., Diez, P. (2018), Numerical Modelling of Multi-Phase Multi-Component Reactive Transport in the Earth`s interior, Geophys. J. Int., 212, 345-388.
It is a (so far matlab) code implementing the general strategies described in the paper above. The approach is truly multiphase in the sense that each dynamic phase is explicitly modelled with an individual set of mass, momentum, energy and chemical mass balance equations coupled via interfacial interaction terms. It is also truly multicomponent in the sense that the compositions of the system and its constituent phases are expressed by a full set of fundamental chemical components (e.g. SiO2, Al2O3, MgO, trace elements and isotopes.) rather than proxies. These chemical components evolve, react with and partition into different phases according to an internally consistent thermodynamic model. We combine concepts from Ensemble Averaging and Classical Irreversible Thermodynamics to obtain sets of macroscopic balance equations that describe the evolution of systems governed by multiphase multicomponent reactive transport (MPMCRT). Equilibrium mineral assemblages, their compositions and physical properties, and closure relations for the balance equations are obtained via a ‘dynamic’ Gibbs free-energy minimization procedure (i.e. minimizations are performed on-the-fly as needed by the simulation).
We are currently using this code to simulate disequilibrium processes during the generation and segregation of melt in the mantle, metasomatism and isotopic fractionation among others. This code is constantly being benchmarked and developed by B. Oliveira and J.C. Afonso at Macquarie University. Please contact us if you are interested in using it for doing research...
On the left: a schematic of the philosophy of the approach... where different "gears" are coupled to make the entire machine/code works
On the right: an example from Oliveira et al. (2018) where a thermal plume develops a second "melt" plume on the top due to the melt percolating through the solid matrix at a faster speed than the actual thermal plume.
Upper Mantle Geoid (UMG3)
This is a filtered geoid model that isolates the lithospheric and upper mantle contributions. It is recommended for studies that target these regions of the planet in order to remove unwanted effects from deep density anomalies. It is presented in:
- Salajeghegh, F. and Afonso, J.C. (2023), The upper mantle geoid for lithospheric structure and dynamics. J. Geophys. Res., https://agupubs.onlinelibrary.wiley.com/doi/abs/10.1029/2023JB026397
LithoRef18: A global reference model of the lithosphere and upper mantle from joint inversion and analysis of multiple data sets
Download model
This model was presented in the paper:
- Afonso, J.C., Salajeghegh, F., Szwillus, W., Ebbing, J. Gaina, C. (2019), A global reference model of the lithosphere and upper mantle from joint inversion and analysis of multiple datasets. Geophys. J. Int., ggz094, https://doi.org/10.1093/gji/ggz094.
LithoRef18 is a new global model for the Earth’s lithosphere and upper mantle obtained through a formal joint inversion of 3D gravity anomalies, geoid height, satellite-derived gravity gradients and absolute elevation complemented with seismic, thermal and petrological prior information . The model includes crustal thickness, average crustal density, lithospheric thickness, depth-dependent density of the lithospheric mantle, lithospheric geotherms, and average density of the sublithospheric mantle down to 410 km depth with a surface discretization of 2°×2°. It is simple, but good...
We also provide some matlab codes to compute far-field effects and interrogate LithoRef18 (see below)
LithoRef18 is also included in the Virtual Globe of the GPLates Portal. Check it out.
"LithoRef model" contains the global model as a text file
"LithoRef code" is a matlab code to compute far-field effects.
This code (cube.m) will allow you to extract profiles along any direction at any location.
Probabilistic inversion of major and trace elements for source composition and physical state
The chemical composition of erupted basalts provides a record of the thermo-chemical state of their source region and the melting conditions that lead to their formation. In Oliveira et al. (2021) we present the first probabilistic inversion framework capable of inverting both trace and major element data of mafic volcanic rocks to constrain mantle potential temperature, depth of melting, and major and trace element source composition. The inversion strategy is based on the combination of i) a two-phase multi-component reactive transport model, ii) a thermodynamic solver for the evolution of major elements and mineral/liquid phases, (iii) a disequilibrium model of trace element partitioning and iv) an adaptive Markov chain Monte Carlo algorithm.
The entire code to perform these type of inversions can be found here:
https://zenodo.org/record/5163125#.YTaTd44zYuV
The associated paper is:
- Oliveira, B., Afonso, J.C. and Kloecking, M. (2021), Melting conditions and mantle source composition from probabilistic joint inversion of major and rare earth element concentrations. Geochim. Cosmochim. Acta (https://eartharxiv.org/repository/view/2606/)
Inversion results for a set of real samples recovered from the Rio Grande region in the US. The panels and histograms depict the posterior probability distribution for relevant model parameters. "ztop" and "Tp" are the bottom of the thermal lithosphere and the mantle potential temperature, respectively.
Some teaching/toy codes
Free-energy minimization using the simplex/Perple_x approach
This simple matlab code implements the main idea of doing Gibbs free energy minimization using linear programming, an approach used by the famous code PerPle_X (e.g. Connolly, 2005). The example allows the user to use different equations of state and mineral phases. The default thermodynamic data is that of Fabrychnaya (1995) for the "olivine" system.
Parallel Tempering, Multiple-try Metropolis sampling
This matlab code illustrates the implementation of the Parallel Tempering (PT) and Mutliple-try Metropolis (MTM) algorithms in the context of MCMC. The user can use PT or MTM separately or together.
Delayed Rejection Adaptive Metropolis (DRAM) in Fortran
At some point I decided to write a simple teaching code in Fortran (but it's very useful for real inversions too!) of the DRAM method of Haario et al. (2006) based on the matlab implementation by Marko Laine.
Approximate Bayesian Computation (ABC)
This is a matlab code that implements ABC to sample from a target distribution. It was written by T. Connell during his MRes project and it is based on the DRAM matlab code by M. Laine.
Gravity modelling and inversion in 2D
A very simple matlab code to illustrate the basics of gravity modelling and inversion using the original algorithm of Talwani et al. (1977). The inversion method is the so-called Compact Gravity Inversion (REF?).
More coming soon...
Make sure you visit this website in a desktop as the mobile version does not contain all software.
I also recommend visiting the websites of the following books, which contain lots of simple matlab examples for Earth sciences:
- Garya, T. (2012) Introduction to Numerical Geodynamic Modelling, CUP.
- Simpson, G. (2016) Practical Finite Element Modelling in earth Sciences using Matlab, Wiley.