diff --git a/contrib/creating-world-notebook/GWB_module.ipynb b/contrib/creating-world-notebook/GWB_module.ipynb new file mode 100644 index 000000000..b1cf31f4e --- /dev/null +++ b/contrib/creating-world-notebook/GWB_module.ipynb @@ -0,0 +1,434 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": { + "id": "LZ7tKD3XvMgY" + }, + "source": [ + "# Creating World with the Geodynamic World Builder \n", + "\n", + "\n", + "\n", + "## Geodynamic World Builder\n", + "This notebook describes the workflow for using the [Geodynamic World Builder](https://geodynamicworldbuilder.github.io/) (GWB) code, a community software used to create complex parameterized geometries to represent the Earth's features such as faults, plumes, and tectonic plates. These features can then be used to described initial conditions, such as temperature or composition, for a geodynamical setting.\n", + "\n", + "### Prerequisites\n", + "- Python 3.1 > and its libraries (json, plotly, meshio)\n", + "- GWB (version 0.5 or later)\n", + "- cmake (version 3.5 or later)\n", + "\n", + "### Installation\n", + "To install GWB, run the following commands in a new cell:\n", + "\n", + "```\n", + "%cd ../../\n", + "! pwd\n", + "! mkdir build\n", + "%cd build\n", + "! pwd\n", + "! cmake .. ; make\n", + "%cd ../contrib/creating-world-notebook/\n", + "\n", + "```\n", + "\n", + "> Note that the above commands changes the current working directory of the Jupyter notebook so make sure that you are in the parent directory of the module after these commands.\n", + "\n", + "More details on the installation can be found in the [manual here](https://gwb.readthedocs.io/en/v1.0.0/user_manual/installation/stand_alone_install.html)" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "id": "MWw_gmO3M2LY" + }, + "source": [ + "## Using GWB\n", + "\n", + "To create your world using the GWB, an input text-based JSON file describing the geometry and properties such as temperature, composition etc. is needed. This file, subsequently referred to as the *World Builder (.wb)* file, is enough to create the complicated geometries with specified initial conditions.\n", + "\n", + "To query the generated world at specific points or visualize it on a specified grid, an additional text file usually with a chosen extension of `.dat` and `.grid` with the executable files `WorldBuilder/build/bin/gwb-dat` and `WorldBuilder/build/bin/gwb-grid` , respectively, is needed.\n", + "\n", + "GWB executables are run from the command line. For example, for a given `my_world.wb` file and a `query-points.dat`, use the command:\n", + "\n", + "`WorldBuilder/build/bin/gwb-dat my_world.wb query-points.dat`\n", + "\n", + "> NOTE : The above line assumes that the executable is in the current working directory. If this is not the case, specify the path to the executable. " + ] + }, + { + "cell_type": "markdown", + "metadata": { + "id": "R4G5a1avUSCO" + }, + "source": [ + "## Writing the input file \n", + "\n", + "All the parameters are described in the GWB manual [here](https://gwb.readthedocs.io/en/v1.0.0/GWB_parameter_listings/world_builder_file/index.html). We briefly describe the parameters by loading the world-builder schema file into a variable, `gwb_schema`.\n", + "\n", + "1. **Global header information** : There are several parameters that describe the global parameters which will be used while creating the World, such as the `version`, `coordinate system` etc. To check the names and description of all of these parameters, you can access the `properties` item stored in the world-builder schema file and the structure for each of the `properties` is defined using `$schema` as `gwb_schema['properties']['$schema']`.\n", + "\n", + "2. **Features** : GWB supports several different feature models corresponding to a [tectonic feature](https://gwb.readthedocs.io/en/v1.0.0/user_manual/parameter_documentation/features/index.html). To check the available features, you can run `gwb_schema['properties']['features']['items']['defaultSnippets']` in this notebook. More information about each of the available features can be found using `gwb_schema['properties']['features']['items']['oneOf'][$num]`, where `$num` is an integer corresponding to the index in the order of the list of models obtained from the previous command. \n", + "\n", + "The geometry of the feature can vary laterally as well as along the down-dip direction for the `fault` and the `subducting plate` models by adding `sections` and `segments`, respectively. *Note that all 'sections' should have the same number of the 'segments'.*\n", + "\n", + "\n", + "Each of these features and their segments and sections can have different properties, such as temperature, composition, grains, velocity. The values of these prescribed properties can be `uniform` or follow some other available distribution for that property. You can get more information on the available property distributions using the schema file. \n", + "\n", + "For example, `gwb_schema['properties']['features']['items']['oneOf'][1]['properties']['composition models']['items']['oneOf']` outputs the list of available composition model for a fault (fault is second in the list of available feature models, therefore we use index 1 to output its properties)." + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": { + "executionInfo": { + "elapsed": 155, + "status": "ok", + "timestamp": 1736877223495, + "user": { + "displayName": "Arushi Saxena", + "userId": "00683773484938428447" + }, + "user_tz": 480 + }, + "id": "jbGqlaVa2jSr" + }, + "outputs": [], + "source": [ + "import subprocess\n", + "import meshio\n", + "from glob import glob\n", + "import plotly.graph_objects as go\n", + "from scipy.interpolate import griddata\n", + "import numpy as np\n", + "\n", + "import json\n", + "from jsonschema import validate" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "# The `.wb` file follows a json schema describing the global parameters and the model features.\n", + "# The schema can be used to investigate the available model features how to describe them \n", + "# in GWB.\n", + "\n", + "schema_file = open('../../doc/world_builder_declarations.schema.json')\n", + "gwb_schema = json.load(schema_file)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Playtime : Creating continents \n", + "\n", + "In this module, we will create continents with their traces following the letters \"GWB\". This example is chosen just to demonstrate how to incorporate a complex geometry using GWB." + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [], + "source": [ + "# STEP 1: First let's create the skeleton file that includes the global header information \n", + "# and the model features. \n", + "\n", + "wb_file = { \"version\": \"1.1\",\n", + " \"interpolation\": \"continuous monotone spline\",\n", + " \"coordinate system\":{\"model\":\"cartesian\"},\n", + " \"features\": []\n", + " }" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [], + "source": [ + "def make_plate_feature(plate_name, coordinates, min_depth=-1, max_depth=10):\n", + " \"\"\"\n", + " This function takes an input fault names and its trace coordinates to create a \n", + " fault feature model for the WB file. \n", + " Several parameters to describe the fault are set to default values for now.\n", + "\n", + " Parameters\n", + " ----------\n", + " fault_name : str\n", + " The name we want to give to the fault.\n", + " coordinates : list\n", + " The coordinates of the fault trace as list of [x, y] points.\n", + " \n", + " Returns\n", + " -------\n", + " out : dict\n", + " The dictionary containing the fault feature model.\n", + " \"\"\"\n", + " \n", + " out = {\n", + " \"model\": \"continental plate\",\n", + " \"name\": plate_name,\n", + " \"min depth\": min_depth,\n", + " \"max depth\": max_depth,\n", + " \"coordinates\": coordinates,\n", + " \n", + " \"composition models\": [{\"model\":\"uniform\", \"compositions\":[0], \"max depth\": 10}] \n", + " }\n", + " \n", + " return out" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [], + "source": [ + "# STEP 2 : Load the accompanying file that contains plate coordinates such that \n", + "# each enclosed plate is separated by a \"<\" character. \n", + "# Store all the feature models as a list of dictionaries, as needed by the GWB\n", + "\n", + "input_file_name = 'plates.csv'\n", + "all_plates_json = []\n", + "\n", + "# Reading the content of the file\n", + "with open(input_file_name, 'r') as file:\n", + " lines = file.read().split(\"<\")\n", + "\n", + "# Determine the number of plates (assuming each block between \"<\" corresponds to one plate)\n", + "number_of_plates = len(lines)\n", + "\n", + "for i in range (1, int(number_of_plates) - 1):\n", + " coordinates = lines[i].strip().split()\n", + "\n", + " # Convert the list of string coordinates to floats\n", + " list_of_floats = [float(item) for item in coordinates]\n", + "\n", + " # Split into x and y coordinates\n", + " coordinates_x = list_of_floats[0::2]\n", + " coordinates_y = list_of_floats[1::2]\n", + "\n", + " # Combine into a coordinate array\n", + " coordinates_array = list(zip(coordinates_x, coordinates_y))\n", + "\n", + " # Append the result\n", + " all_plates_json.append(make_plate_feature('plate_' + str(i), coordinates_array))" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [], + "source": [ + "# STEP 3: Create the json file that includes all the continental plates \n", + "\n", + "output_wb_file = \"my_world_builder.wb\"\n", + "wb_file['features'] = all_plates_json\n", + "\n", + "with open(output_wb_file, \"w\") as f :\n", + " json.dump(wb_file, f)" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + " \r" + ] + }, + { + "data": { + "text/plain": [ + "CompletedProcess(args='../../build/bin/gwb-grid my_world_builder.wb cartesian_surface.grid', returncode=0)" + ] + }, + "execution_count": 7, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "# STEP 4: Run the GWB using the .wb file generated above using the\n", + "# accompanying .grid file, 'cartesian_surface.grid'\n", + "# This .grid file is a simple text file that contains the extent and\n", + "# the discretization of the grid. We need to make sure that the extent \n", + "# of the grid is big enough to include all the plates.\n", + "\n", + "gwb_grid_location = \"../../build/bin/gwb-grid \"\n", + "grid_file_name = \" cartesian_surface.grid\"\n", + "command_to_run = gwb_grid_location + output_wb_file + grid_file_name\n", + "\n", + "subprocess.run(command_to_run, shell=True)" + ] + }, + { + "cell_type": "markdown", + "metadata": { + "id": "iaW38MIj2a_i" + }, + "source": [ + "## Visualizing the results\n", + "\n", + "### Meshio\n", + "To visualize the generated world, we use the `gwb-grid` to query the properties on a grid file, `cartesian_surface.grid` available as one of the [tests](https://github.com/GeodynamicWorldBuilder/WorldBuilder/tree/main/tests) in GWB." + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": { + "id": "RqKJjV7K6w8H" + }, + "outputs": [], + "source": [ + "def read_vtu_file (file_name):\n", + " \"\"\"\n", + " This function takes an input .vtu file_name (either structured or\n", + " unstructured mesh) and returns the point coordinates and values in that file. \n", + " \n", + " Parameters\n", + " ----------\n", + " file_name : str\n", + " The name of the vtu file\n", + " \n", + " Returns\n", + " -------\n", + " nodes : A numpy array with columns corresponding to the \n", + " coordinates of the node points (x, y, z) in the input mesh.\n", + " \n", + " data_out : A dictionary containing key:values corresponding to the \n", + " available fields and their values in the input .vtu file.\n", + " \"\"\"\n", + " \n", + " mesh = meshio.read(file_name)\n", + " nodes = mesh.points\n", + " data_out = mesh.point_data\n", + " \n", + " return (nodes, data_out)" + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": { + "executionInfo": { + "elapsed": 145, + "status": "ok", + "timestamp": 1736876702953, + "user": { + "displayName": "Arushi Saxena", + "userId": "00683773484938428447" + }, + "user_tz": 480 + }, + "id": "_QSpq5hL6YrG" + }, + "outputs": [], + "source": [ + "vtu_file = output_wb_file.replace('.wb', '.vtu')\n", + "coordinates, data = read_vtu_file (vtu_file)" + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "metadata": {}, + "outputs": [], + "source": [ + "%%capture \n", + "\n", + "# Create the 3D scatter plot using Plotly\n", + "points_x = coordinates[:, 0]\n", + "points_y = coordinates[:, 1]\n", + "points_z = coordinates[:, 2]\n", + "\n", + "x_ = np.linspace(min(points_x), max(points_x), 80)\n", + "y_ = np.linspace(min(points_y), max(points_y), 40)\n", + "z_ = np.linspace(min(points_z), max(points_z), 5)\n", + "\n", + "values = data['Composition 0']\n", + "\n", + "X, Y, Z = np.meshgrid(x_, y_, z_, indexing='ij')\n", + "values = griddata((points_x, points_y, points_z), values, (X, Y, Z), method='nearest')\n", + "\n", + "# Create a meshgrid with the spherical coordinates\n", + "\n", + "fig = go.Figure(data=go.Volume(\n", + " x=X.flatten(),\n", + " y=Y.flatten(),\n", + " z=Z.flatten(),\n", + " value=values.flatten(), \n", + " colorscale='rdbu_r', \n", + " colorbar = dict(title='Composition', tickvals=[0, 1], ticktext=['0', '1']),\n", + " opacity=0.2,\n", + " surface_count=30,\n", + " ))\n", + "\n", + "camera = dict( eye=dict(x=0, y=0, z=4), \n", + " up=dict(x=0, y=1, z=0), \n", + " center=dict(x=0, y=0, z=0))\n", + "\n", + "fig.update_scenes(aspectmode='data')\n", + "fig.update_layout(scene_camera=camera, title=\"GWB Output\")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Paraview Glance\n", + "Another tool to visualize the models is to use the [Glance](https://kitware.github.io/glance/app/) app. You can easily load the generated `.vtu` file from the GWB in Glance and rotate around the modeled results and modify the visualization field." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Exercise: things to try\n", + "\n", + "* Change the max depth \n", + "* Change composition value to 1 and modify the current .grid file for 2 compositions\n", + "* Plot temperatures instead of Composition\n", + "* Add more features\n", + "* Modify the resolution by changing, one of more of n_cell_x, n_cell_y, n_cell_z, in the.grid file" + ] + } + ], + "metadata": { + "colab": { + "authorship_tag": "ABX9TyPTW1UuXYKrsm8pvCCZ0sX7", + "provenance": [] + }, + "kernelspec": { + "display_name": "base", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.11.5" + } + }, + "nbformat": 4, + "nbformat_minor": 1 +} diff --git a/contrib/creating-world-notebook/cartesian_surface.grid b/contrib/creating-world-notebook/cartesian_surface.grid new file mode 100644 index 000000000..1d76fd7c9 --- /dev/null +++ b/contrib/creating-world-notebook/cartesian_surface.grid @@ -0,0 +1,19 @@ +# output variables +grid_type = cartesian +dim = 3 +compositions = 1 #5 +vtu_output_format = ASCII + +# domain of the grid +x_min = -10 +x_max = 1250 +y_min = 0 +y_max = 410 +z_min = -1 +z_max = 10 + + +# grid properties +n_cell_x = 120 +n_cell_y = 40 +n_cell_z = 10 diff --git a/contrib/creating-world-notebook/plates.csv b/contrib/creating-world-notebook/plates.csv new file mode 100644 index 000000000..172c5262a --- /dev/null +++ b/contrib/creating-world-notebook/plates.csv @@ -0,0 +1,254 @@ +< +324.899 375.768 +314.165 381.317 +314.165 381.317 +300.118 387.123 +285.237 392.206 +267.037 396.760 +247.998 399.806 +233.940 404.039 +215.712 404.662 +195.017 407.836 +171.808 407.271 +146.938 406.049 +125.358 401.425 +103.778 396.802 +84.666 389.629 +68.851 379.841 +53.024 368.480 +35.491 350.174 +25.446 336.007 +25.446 336.007 +18.697 319.225 +8.629 301.914 +5.154 279.374 +3.358 259.850 +-0.095 240.454 +3.079 220.546 +2.116 201.744 +5.295 182.621 +15.873 155.062 +24.828 132.347 +37.157 118.022 +51.115 99.639 +74.195 82.124 +107.225 65.412 +134.519 57.796 +170.913 47.902 +204.883 46.848 +248.804 46.597 +278.650 48.221 +308.507 51.418 +327.591 54.661 +327.980 226.092 +260.901 232.853 +207.069 237.803 +191.326 238.234 +192.709 199.587 +239.915 195.150 +275.536 193.181 +276.158 164.032 +279.248 132.333 +279.041 103.247 +276.450 88.504 +326.482 131.825 +326.280 103.526 +327.673 182.857 +288.234 230.739 +246.582 83.735 +206.814 85.238 +159.619 91.248 +126.555 103.243 +104.264 115.192 +80.395 138.274 +67.327 165.240 +56.749 192.799 +55.254 215.724 +54.627 244.087 +53.188 274.872 +67.459 300.510 +81.713 323.790 +100.902 341.969 +131.666 356.106 +159.071 364.211 +196.383 366.831 +236.118 360.612 +276.642 348.826 +303.925 339.637 +314.848 360.815 +< +374.702 406.860 +382.050 391.348 +386.892 372.884 +390.110 359.264 +396.619 342.244 +400.643 325.416 +407.941 302.829 +415.200 274.739 +423.282 245.800 +429.707 216.988 +441.096 187.006 +445.893 162.253 +456.437 129.977 +462.923 109.812 +471.000 80.087 +478.264 52.783 +497.303 49.738 +521.318 47.094 +530.423 45.603 +536.355 64.021 +544.844 92.465 +540.614 80.208 +549.942 110.161 +559.221 133.039 +566.065 163.184 +565.981 151.392 +575.288 178.201 +580.369 193.538 +586.290 210.383 +593.056 229.523 +598.204 254.293 +607.427 269.310 +612.536 288.578 +619.324 310.862 +624.444 331.702 +634.272 315.212 +637.423 292.159 +642.253 272.123 +648.751 253.530 +654.404 232.643 +660.885 211.693 +668.989 185.897 +677.910 158.466 +689.349 135.559 +692.534 117.223 +700.633 90.641 +708.765 68.776 +719.325 38.858 +725.073 31.335 +751.589 30.857 +768.988 30.298 +776.582 49.374 +783.393 74.802 +792.699 101.610 +804.535 134.515 +798.609 116.884 +811.340 159.158 +821.469 185.116 +826.590 205.956 +833.383 229.026 +841.022 254.390 +848.639 276.610 +857.973 307.349 +865.628 335.071 +873.245 357.291 +877.464 367.976 +849.298 369.368 +831.904 370.713 +821.959 370.696 +816.805 345.139 +810.000 320.497 +803.234 301.358 +798.925 278.095 +790.458 252.795 +783.687 232.869 +774.358 202.917 +765.891 177.616 +761.594 155.926 +754.811 134.428 +750.569 120.599 +746.310 104.412 +745.376 89.540 +733.209 126.661 +725.955 155.537 +717.879 185.262 +706.512 218.389 +694.334 253.937 +684.601 283.791 +676.525 313.517 +664.313 344.349 +655.398 372.567 +651.346 385.464 +625.686 389.808 +604.979 391.410 +598.224 373.842 +591.486 358.633 +586.388 340.938 +580.484 326.450 +573.741 310.455 +567.809 292.037 +560.221 273.748 +553.466 256.180 +548.323 232.196 +539.028 206.960 +526.447 185.911 +524.607 160.098 +517.007 140.236 +504.443 121.545 +500.508 150.950 +493.216 174.323 +481.843 206.664 +470.549 250.009 +460.839 283.007 +451.246 332.513 +439.935 373.500 +429.379 404.204 +397.088 408.275 +475.396 232.331 +< +938.762 389.022 +939.484 374.022 +939.373 358.300 +939.283 345.723 +937.471 323.841 +937.297 299.472 +937.992 280.542 +936.996 257.024 +936.834 234.227 +937.495 210.580 +937.310 184.640 +937.160 163.415 +937.009 142.191 +936.892 125.683 +937.608 109.898 +938.291 89.395 +935.673 70.721 +936.400 56.508 +937.123 41.508 +962.771 35.592 +988.453 34.392 +1017.443 32.150 +1045.621 32.330 +1081.237 29.576 +1110.226 27.334 +1135.942 30.850 +1160.863 39.148 +1182.498 51.632 +1195.012 63.248 +1206.736 80.431 +1212.724 106.710 +1211.246 131.993 +1203.970 157.724 +1186.727 180.293 +1175.193 189.837 +1155.376 200.021 +1140.483 203.533 +1128.087 208.424 +1149.695 216.978 +1167.984 225.001 +1186.301 236.956 +1194.695 252.037 +1200.661 275.171 +1200.028 302.748 +1188.639 332.730 +1172.179 348.947 +1152.368 359.917 +1124.268 370.742 +1105.229 373.788 +1083.699 376.239 +1053.887 379.331 +1024.080 383.210 +994.262 385.516 +973.555 387.117 +954.505 388.591 +<