diff --git a/examples/symba_mars_disk/param.in b/examples/symba_mars_disk/param.in index a50be00e4..725927f1b 100644 --- a/examples/symba_mars_disk/param.in +++ b/examples/symba_mars_disk/param.in @@ -1,13 +1,13 @@ !Parameter file for the SyMBA-RINGMOONS test T0 0.0 -TSTOP 6000.0 +TSTOP 60000.0 DT 600.0 CB_IN cb.in PL_IN mars.in TP_IN tp.in IN_TYPE ASCII -ISTEP_OUT 1 -ISTEP_DUMP 1 +ISTEP_OUT 10 +ISTEP_DUMP 10 !BIN_OUT bin.dat !OUT_TYPE REAL8 BIN_OUT bin.nc diff --git a/examples/symba_mars_disk/testnetcdf.ipynb b/examples/symba_mars_disk/testnetcdf.ipynb index 955dc994f..a6d21ec0e 100644 --- a/examples/symba_mars_disk/testnetcdf.ipynb +++ b/examples/symba_mars_disk/testnetcdf.ipynb @@ -35,9 +35,30 @@ "text": [ "Reading Swiftest file param.in\n", "\n", - "Creating Dataset\n", - "Successfully converted 11 output frames.\n", - "Swiftest simulation data stored as xarray DataSet .ds\n" + "Creating Dataset\n" + ] + }, + { + "ename": "UnicodeDecodeError", + "evalue": "'utf-8' codec can't decode byte 0x80 in position 0: invalid start byte", + "output_type": "error", + "traceback": [ + "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", + "\u001b[0;31mUnicodeDecodeError\u001b[0m Traceback (most recent call last)", + "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[0msim\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mswiftest\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mSimulation\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mparam_file\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m\"param.in\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 2\u001b[0;31m \u001b[0msim\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mbin2xr\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m", + "\u001b[0;32m~/git/swiftest/python/swiftest/swiftest/simulation_class.py\u001b[0m in \u001b[0;36mbin2xr\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m 172\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0mbin2xr\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 173\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mcodename\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0;34m\"Swiftest\"\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 174\u001b[0;31m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mds\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mio\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mswiftest2xr\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mparam\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 175\u001b[0m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m'Swiftest simulation data stored as xarray DataSet .ds'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 176\u001b[0m \u001b[0;32melif\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mcodename\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0;34m\"Swifter\"\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/git/swiftest/python/swiftest/swiftest/io.py\u001b[0m in \u001b[0;36mswiftest2xr\u001b[0;34m(param)\u001b[0m\n\u001b[1;32m 697\u001b[0m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m'\\nCreating Dataset'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 698\u001b[0m \u001b[0mds\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mxr\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mopen_dataset\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mparam\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m'BIN_OUT'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 699\u001b[0;31m \u001b[0mds\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mclean_string_values\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mparam\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mds\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 700\u001b[0m \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 701\u001b[0m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34mf\"Error encountered. OUT_TYPE {param['OUT_TYPE']} not recognized.\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/git/swiftest/python/swiftest/swiftest/io.py\u001b[0m in \u001b[0;36mclean_string_values\u001b[0;34m(param, ds)\u001b[0m\n\u001b[1;32m 719\u001b[0m \u001b[0mds\u001b[0m \u001b[0;34m:\u001b[0m \u001b[0mxarray\u001b[0m \u001b[0mdataset\u001b[0m \u001b[0;32mwith\u001b[0m \u001b[0mthe\u001b[0m \u001b[0mstrings\u001b[0m \u001b[0mcleaned\u001b[0m \u001b[0mup\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 720\u001b[0m \"\"\" \n\u001b[0;32m--> 721\u001b[0;31m \u001b[0mds\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m'name'\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mds\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m'name'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mstr\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdecode\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mencoding\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m'utf-8'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 722\u001b[0m \u001b[0mds\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m'particle_type'\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mds\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m'particle_type'\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mstr\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdecode\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mencoding\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m'utf-8'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 723\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0;34m'origin_type'\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mds\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/.conda/envs/cent7/2020.02-py37/swiftestOOF/lib/python3.7/site-packages/xarray/core/accessor_str.py\u001b[0m in \u001b[0;36mdecode\u001b[0;34m(self, encoding, errors)\u001b[0m\n\u001b[1;32m 2545\u001b[0m \u001b[0mdecoder\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mcodecs\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgetdecoder\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mencoding\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2546\u001b[0m \u001b[0mfunc\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;32mlambda\u001b[0m \u001b[0mx\u001b[0m\u001b[0;34m:\u001b[0m \u001b[0mdecoder\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mx\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0merrors\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 2547\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_apply\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mfunc\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mfunc\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mdtype\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mstr_\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 2548\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2549\u001b[0m def encode(\n", + "\u001b[0;32m~/.conda/envs/cent7/2020.02-py37/swiftestOOF/lib/python3.7/site-packages/xarray/core/accessor_str.py\u001b[0m in \u001b[0;36m_apply\u001b[0;34m(self, func, dtype, output_core_dims, output_sizes, func_args, func_kwargs)\u001b[0m\n\u001b[1;32m 239\u001b[0m \u001b[0moutput_sizes\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0moutput_sizes\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 240\u001b[0m \u001b[0mfunc_args\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mfunc_args\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 241\u001b[0;31m \u001b[0mfunc_kwargs\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mfunc_kwargs\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 242\u001b[0m )\n\u001b[1;32m 243\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/.conda/envs/cent7/2020.02-py37/swiftestOOF/lib/python3.7/site-packages/xarray/core/accessor_str.py\u001b[0m in \u001b[0;36m_apply_str_ufunc\u001b[0;34m(func, obj, dtype, output_core_dims, output_sizes, func_args, func_kwargs)\u001b[0m\n\u001b[1;32m 136\u001b[0m \u001b[0moutput_core_dims\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0moutput_core_dims\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 137\u001b[0m \u001b[0mdask_gufunc_kwargs\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mdask_gufunc_kwargs\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 138\u001b[0;31m \u001b[0;34m**\u001b[0m\u001b[0mfunc_kwargs\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 139\u001b[0m )\n\u001b[1;32m 140\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/.conda/envs/cent7/2020.02-py37/swiftestOOF/lib/python3.7/site-packages/xarray/core/computation.py\u001b[0m in \u001b[0;36mapply_ufunc\u001b[0;34m(func, input_core_dims, output_core_dims, exclude_dims, vectorize, join, dataset_join, dataset_fill_value, keep_attrs, kwargs, dask, output_dtypes, output_sizes, meta, dask_gufunc_kwargs, *args)\u001b[0m\n\u001b[1;32m 1142\u001b[0m \u001b[0mjoin\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mjoin\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 1143\u001b[0m \u001b[0mexclude_dims\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mexclude_dims\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 1144\u001b[0;31m \u001b[0mkeep_attrs\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mkeep_attrs\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 1145\u001b[0m )\n\u001b[1;32m 1146\u001b[0m \u001b[0;31m# feed Variables directly through apply_variable_ufunc\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/.conda/envs/cent7/2020.02-py37/swiftestOOF/lib/python3.7/site-packages/xarray/core/computation.py\u001b[0m in \u001b[0;36mapply_dataarray_vfunc\u001b[0;34m(func, signature, join, exclude_dims, keep_attrs, *args)\u001b[0m\n\u001b[1;32m 269\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 270\u001b[0m \u001b[0mdata_vars\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0mgetattr\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0ma\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m\"variable\"\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0ma\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0ma\u001b[0m \u001b[0;32min\u001b[0m \u001b[0margs\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 271\u001b[0;31m \u001b[0mresult_var\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfunc\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0mdata_vars\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 272\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 273\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0msignature\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mnum_outputs\u001b[0m \u001b[0;34m>\u001b[0m \u001b[0;36m1\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/.conda/envs/cent7/2020.02-py37/swiftestOOF/lib/python3.7/site-packages/xarray/core/computation.py\u001b[0m in \u001b[0;36mapply_variable_ufunc\u001b[0;34m(func, signature, exclude_dims, dask, output_dtypes, vectorize, keep_attrs, dask_gufunc_kwargs, *args)\u001b[0m\n\u001b[1;32m 722\u001b[0m )\n\u001b[1;32m 723\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 724\u001b[0;31m \u001b[0mresult_data\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfunc\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0minput_data\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 725\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 726\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0msignature\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mnum_outputs\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0;36m1\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/.local/lib/python3.7/site-packages/numpy/lib/function_base.py\u001b[0m in \u001b[0;36m__call__\u001b[0;34m(self, *args, **kwargs)\u001b[0m\n\u001b[1;32m 2111\u001b[0m \u001b[0mvargs\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mextend\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mkwargs\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0m_n\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0m_n\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mnames\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2112\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 2113\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_vectorize_call\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mfunc\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mfunc\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0margs\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mvargs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 2114\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2115\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0m_get_ufunc_and_otypes\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mfunc\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0margs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/.local/lib/python3.7/site-packages/numpy/lib/function_base.py\u001b[0m in \u001b[0;36m_vectorize_call\u001b[0;34m(self, func, args)\u001b[0m\n\u001b[1;32m 2195\u001b[0m for a in args]\n\u001b[1;32m 2196\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 2197\u001b[0;31m \u001b[0moutputs\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mufunc\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0minputs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 2198\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2199\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mufunc\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mnout\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0;36m1\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/.conda/envs/cent7/2020.02-py37/swiftestOOF/lib/python3.7/site-packages/xarray/core/accessor_str.py\u001b[0m in \u001b[0;36m\u001b[0;34m(x)\u001b[0m\n\u001b[1;32m 2541\u001b[0m \"\"\"\n\u001b[1;32m 2542\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mencoding\u001b[0m \u001b[0;32min\u001b[0m \u001b[0m_cpython_optimized_decoders\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 2543\u001b[0;31m \u001b[0mfunc\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;32mlambda\u001b[0m \u001b[0mx\u001b[0m\u001b[0;34m:\u001b[0m \u001b[0mx\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdecode\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mencoding\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0merrors\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 2544\u001b[0m \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2545\u001b[0m \u001b[0mdecoder\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mcodecs\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mgetdecoder\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mencoding\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;31mUnicodeDecodeError\u001b[0m: 'utf-8' codec can't decode byte 0x80 in position 0: invalid start byte" ] } ], @@ -50,388 +71,14 @@ "cell_type": "code", "execution_count": 3, "metadata": {}, - "outputs": [ - { - "data": { - "text/html": [ - "
\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "\n", - "
<xarray.DataArray 'origin_type' (id: 1521)>\n",
-       "array(['Initial conditions', 'Body                    Initial',\n",
-       "       'Body                    Initial', ..., '', '', ''], dtype='<U32')\n",
-       "Coordinates:\n",
-       "  * id       (id) int32 0 1 2 3 4 5 6 7 ... 1514 1515 1516 1517 1518 1519 1520
" - ], - "text/plain": [ - "\n", - "array(['Initial conditions', 'Body Initial',\n", - " 'Body Initial', ..., '', '', ''], dtype='
<xarray.Dataset>\n",
-       "Dimensions:        ()\n",
+       "
<xarray.DataArray 'name' ()>\n",
+       "array(b'UNNAMED', dtype='|S32')\n",
        "Coordinates:\n",
-       "    time           float64 0.0\n",
-       "    id             int32 2\n",
-       "Data variables: (12/48)\n",
-       "    npl            int32 1500\n",
-       "    ntp            int32 0\n",
-       "    name           <U32 'Body3'\n",
-       "    particle_type  <U12 'Massive Body'\n",
-       "    xhx            float64 -8.012e+06\n",
-       "    xhy            float64 -6.936e+06\n",
-       "    ...             ...\n",
-       "    origin_xhx     float64 6.013e-154\n",
-       "    origin_xhy     float64 6.013e-154\n",
-       "    origin_xhz     float64 0.0\n",
-       "    origin_vhx     float64 7.029e+06\n",
-       "    origin_vhy     float64 6.052e+06\n",
-       "    origin_vhz     float64 3.879e+03
"
+       "    id       int32 1885
" ], "text/plain": [ - "\n", - "Dimensions: ()\n", + "\n", + "array(b'UNNAMED', dtype='|S32')\n", "Coordinates:\n", - " time float64 0.0\n", - " id int32 2\n", - "Data variables: (12/48)\n", - " npl int32 ...\n", - " ntp int32 ...\n", - " name netcdf_close !! Closes an open NetCDF file procedure :: initialize => netcdf_initialize_output !! Initialize a set of parameters used to identify a NetCDF output object @@ -136,7 +145,7 @@ module swiftest_classes logical :: lyarkovsky = .false. !! Turn on Yarkovsky effect logical :: lyorp = .false. !! Turn on YORP effect - class(netcdf_parameters), allocatable :: nciu !! Object containing NetCDF parameters + type(netcdf_parameters) :: nciu !! Object containing NetCDF parameters contains procedure :: reader => io_param_reader procedure :: writer => io_param_writer @@ -153,6 +162,10 @@ module swiftest_classes type :: swiftest_particle_info character(len=NAMELEN) :: name !! Non-unique name character(len=NAMELEN) :: particle_type !! String containing a description of the particle type (e.g. Central Body, Massive Body, Test Particle) + character(len=NAMELEN) :: origin_type !! String containing a description of the origin of the particle (e.g. Initial Conditions, Supercatastrophic, Disruption, etc.) + real(DP) :: origin_time !! The time of the particle's formation + real(DP), dimension(NDIM) :: origin_xh !! The heliocentric distance vector at the time of the particle's formation + real(DP), dimension(NDIM) :: origin_vh !! The heliocentric velocity vector at the time of the particle's formation contains procedure :: dump => io_dump_particle_info !! Dumps contents of particle information to file procedure :: read_in => io_read_in_particle_info !! Read in a particle information object from an open file @@ -176,7 +189,7 @@ module swiftest_classes !******************************************************************************************************************************** !> A concrete lass for the central body in a Swiftest simulation type, abstract, extends(swiftest_base) :: swiftest_cb - class(swiftest_particle_info), allocatable :: info !! Particle metadata information + type(swiftest_particle_info) :: info !! Particle metadata information integer(I4B) :: id = 0 !! External identifier (unique) real(DP) :: mass = 0.0_DP !! Central body mass (units MU) real(DP) :: Gmass = 0.0_DP !! Central mass gravitational term G * mass (units GU * MU) @@ -215,7 +228,7 @@ module swiftest_classes !! Superclass that defines the generic elements of a Swiftest particle logical :: lfirst = .true. !! Run the current step as a first integer(I4B) :: nbody = 0 !! Number of bodies - class(swiftest_particle_info), dimension(:), allocatable :: info !! Particle metadata information + type(swiftest_particle_info), dimension(:), allocatable :: info !! Particle metadata information integer(I4B), dimension(:), allocatable :: id !! External identifier (unique) integer(I4B), dimension(:), allocatable :: status !! An integrator-specific status indicator logical, dimension(:), allocatable :: ldiscard !! Body should be discarded @@ -1002,7 +1015,7 @@ end subroutine setup_encounter module subroutine setup_initialize_particle_info_system(self, param) implicit none class(swiftest_nbody_system), intent(inout) :: self !! Swiftest nbody system object - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters + class(swiftest_parameters), intent(inout) :: param !! Current run configuration parameters end subroutine setup_initialize_particle_info_system module subroutine setup_initialize_system(self, param) @@ -1084,8 +1097,8 @@ end subroutine util_append_arr_I4B module subroutine util_append_arr_info(arr, source, nold, nsrc, lsource_mask) implicit none - class(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Destination array - class(swiftest_particle_info), dimension(:), allocatable, intent(in) :: source !! Array to append + type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Destination array + type(swiftest_particle_info), dimension(:), allocatable, intent(in) :: source !! Array to append integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to end subroutine util_append_arr_info @@ -1253,8 +1266,8 @@ end subroutine util_fill_arr_I4B module subroutine util_fill_arr_info(keeps, inserts, lfill_list) implicit none - class(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - class(swiftest_particle_info), dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep + type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep + type(swiftest_particle_info), dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps end subroutine util_fill_arr_info @@ -1322,7 +1335,7 @@ end subroutine util_resize_arr_I4B module subroutine util_resize_arr_info(arr, nnew) implicit none - class(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Array to resize + type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Array to resize integer(I4B), intent(in) :: nnew !! New size end subroutine util_resize_arr_info @@ -1506,7 +1519,7 @@ end subroutine util_sort_rearrange_arr_I4B module subroutine util_sort_rearrange_arr_info(arr, ind, n) implicit none - class(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Destination array + type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Destination array integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange end subroutine util_sort_rearrange_arr_info @@ -1603,8 +1616,8 @@ end subroutine util_spill_arr_I8B module subroutine util_spill_arr_info(keeps, discards, lspill_list, ldestructive) implicit none - class(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - class(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: discards !! Array of discards + type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep + type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: discards !! Array of discards logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discardss logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not end subroutine util_spill_arr_info diff --git a/src/modules/swiftest_globals.f90 b/src/modules/swiftest_globals.f90 index b84e7f2f9..a9dd2d762 100644 --- a/src/modules/swiftest_globals.f90 +++ b/src/modules/swiftest_globals.f90 @@ -181,4 +181,14 @@ module swiftest_globals character(*), parameter :: ECOLLISIONS_VARNAME = "Ecollisions" !! NetCDF name of the escaped angular momentum y variable character(*), parameter :: EUNTRACKED_VARNAME = "Euntracked" !! NetCDF name of the energy that is untracked due to loss (untracked potential energy due to mergers and body energy for escaped bodies) character(*), parameter :: GMESCAPE_VARNAME = "GMescape" !! NetCDF name of the G*Mass of bodies that escape the system + character(*), parameter :: ORIGIN_TYPE_VARNAME = "origin_type" + character(*), parameter :: ORIGIN_TIME_VARNAME = "origin_time" + character(*), parameter :: ORIGIN_XHX_VARNAME = "origin_xhx" + character(*), parameter :: ORIGIN_XHY_VARNAME = "origin_xhy" + character(*), parameter :: ORIGIN_XHZ_VARNAME = "origin_xhz" + character(*), parameter :: ORIGIN_VHX_VARNAME = "origin_vhx" + character(*), parameter :: ORIGIN_VHY_VARNAME = "origin_vhy" + character(*), parameter :: ORIGIN_VHZ_VARNAME = "origin_vhz" + character(*), parameter :: PL_TINY_TYPE_NAME = "Semi-Interacting Massive Body" + end module swiftest_globals diff --git a/src/modules/symba_classes.f90 b/src/modules/symba_classes.f90 index 8bc32ffed..3a9e37fce 100644 --- a/src/modules/symba_classes.f90 +++ b/src/modules/symba_classes.f90 @@ -14,29 +14,6 @@ module symba_classes integer(I4B), private, parameter :: NTENC = 3 real(DP), private, parameter :: RHSCALE = 6.5_DP real(DP), private, parameter :: RSHELL = 0.48075_DP - character(*), parameter :: ORIGIN_TYPE_VARNAME = "origin_type" - character(*), parameter :: ORIGIN_TIME_VARNAME = "origin_time" - character(*), parameter :: ORIGIN_XHX_VARNAME = "origin_xhx" - character(*), parameter :: ORIGIN_XHY_VARNAME = "origin_xhy" - character(*), parameter :: ORIGIN_XHZ_VARNAME = "origin_xhz" - character(*), parameter :: ORIGIN_VHX_VARNAME = "origin_vhx" - character(*), parameter :: ORIGIN_VHY_VARNAME = "origin_vhy" - character(*), parameter :: ORIGIN_VHZ_VARNAME = "origin_vhz" - character(*), parameter :: PL_TINY_TYPE_NAME = "Semi-Interacting Massive Body" - - type, extends(netcdf_parameters) :: symba_netcdf_parameters - integer(I4B) :: origin_type_varid !! NetCDF ID for the origin type - integer(I4B) :: origin_time_varid !! NetCDF ID for the origin type - integer(I4B) :: origin_xhx_varid !! NetCDF ID for the origin xh x component - integer(I4B) :: origin_xhy_varid !! NetCDF ID for the origin xh y component - integer(I4B) :: origin_xhz_varid !! NetCDF ID for the origin xh z component - integer(I4B) :: origin_vhx_varid !! NetCDF ID for the origin xh x component - integer(I4B) :: origin_vhy_varid !! NetCDF ID for the origin xh y component - integer(I4B) :: origin_vhz_varid !! NetCDF ID for the origin xh z component - contains - procedure :: initialize => symba_netcdf_initialize_output !! Initialize a set of parameters used to identify a NetCDF output objec - procedure :: open => symba_netcdf_open !! Opens a NetCDF file - end type symba_netcdf_parameters type, extends(swiftest_parameters) :: symba_parameters real(DP) :: GMTINY = -1.0_DP !! Smallest G*mass that is fully gravitating @@ -48,20 +25,6 @@ module symba_classes procedure :: writer => symba_io_param_writer end type symba_parameters - !******************************************************************************************************************************** - ! symba_swiftest_particle_info class definitions and method interfaces - !******************************************************************************************************************************* - !> Class definition for the particle origin information object. This object is used to track time, location, and collisional regime - !> of fragments produced in collisional events. - type, extends(swiftest_particle_info) :: symba_particle_info - character(len=NAMELEN) :: origin_type !! String containing a description of the origin of the particle (e.g. Initial Conditions, Supercatastrophic, Disruption, etc.) - real(DP) :: origin_time !! The time of the particle's formation - real(DP), dimension(NDIM) :: origin_xh !! The heliocentric distance vector at the time of the particle's formation - real(DP), dimension(NDIM) :: origin_vh !! The heliocentric velocity vector at the time of the particle's formation - contains - procedure :: read_in => symba_io_read_in_particle_info !! Reads in SyMBA particle information metadata from an open unformatted file - end type symba_particle_info - !******************************************************************************************************************************** ! symba_kinship class definitions and method interfaces !******************************************************************************************************************************* @@ -82,7 +45,6 @@ module symba_classes real(DP) :: R0 = 0.0_DP !! Initial radius of the central body real(DP) :: dR = 0.0_DP !! Change in the radius of the central body contains - procedure :: write_frame_netcdf => symba_netcdf_write_frame_cb !! I/O routine for writing out a single frame of time-series data for all bodies in the system in NetCDF format end type symba_cb !******************************************************************************************************************************** @@ -120,7 +82,6 @@ module symba_classes procedure :: sort => symba_util_sort_pl !! Sorts body arrays by a sortable componen procedure :: rearrange => symba_util_sort_rearrange_pl !! Rearranges the order of array elements of body based on an input index array. Used in sorting methods procedure :: spill => symba_util_spill_pl !! "Spills" bodies from one object to another depending on the results of a mask (uses the PACK intrinsic) - procedure :: write_frame_netcdf => symba_netcdf_write_frame_pl !! I/O routine for writing out a single frame of time-series data for all bodies in the system in NetCDF format end type symba_pl type, extends(symba_pl) :: symba_merger @@ -150,7 +111,6 @@ module symba_classes procedure :: sort => symba_util_sort_tp !! Sorts body arrays by a sortable componen procedure :: rearrange => symba_util_sort_rearrange_tp !! Rearranges the order of array elements of body based on an input index array. Used in sorting methods procedure :: spill => symba_util_spill_tp !! "Spills" bodies from one object to another depending on the results of a mask (uses the PACK intrinsic) - procedure :: write_frame_netcdf => symba_netcdf_write_frame_tp !! I/O routine for writing out a single frame of time-series data for all bodies in the system in NetCDF format end type symba_tp !******************************************************************************************************************************** @@ -202,7 +162,6 @@ module symba_classes contains procedure :: write_discard => symba_io_write_discard !! Write out information about discarded and merged planets and test particles in SyMBA procedure :: initialize => symba_setup_initialize_system !! Performs SyMBA-specific initilization steps - procedure :: init_particle_info => symba_setup_initialize_particle_info_system !! Initialize the system from input files procedure :: step => symba_step_system !! Advance the SyMBA nbody system forward in time by one step procedure :: interp => symba_step_interp_system !! Perform an interpolation step on the SymBA nbody system procedure :: set_recur_levels => symba_step_set_recur_levels_system !! Sets recursion levels of bodies and encounter lists to the current system level @@ -384,13 +343,6 @@ module subroutine symba_util_index_eucl_plpl(self, param) class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters end subroutine symba_util_index_eucl_plpl - module subroutine symba_io_dump_particle_info(self, iu) - implicit none - class(symba_particle_info), intent(in) :: self !! Particle metadata information object - integer(I4B), intent(in) :: iu !! Open file unit number - end subroutine symba_io_dump_particle_info - - module subroutine symba_io_param_reader(self, unit, iotype, v_list, iostat, iomsg) implicit none class(symba_parameters), intent(inout) :: self !! Current run configuration parameters with SyMBA additionss @@ -413,12 +365,6 @@ module subroutine symba_io_param_writer(self, unit, iotype, v_list, iostat, ioms character(len=*), intent(inout) :: iomsg !! Message to pass if iostat /= 0 end subroutine symba_io_param_writer - module subroutine symba_io_read_in_particle_info(self, iu) - implicit none - class(symba_particle_info), intent(inout) :: self !! Particle metadata information object - integer(I4B), intent(in) :: iu !! Open file unit number - end subroutine symba_io_read_in_particle_info - module subroutine symba_io_write_discard(self, param) use swiftest_classes, only : swiftest_parameters implicit none @@ -460,51 +406,6 @@ module subroutine symba_kick_encounter(self, system, dt, irec, sgn) integer(I4B), intent(in) :: sgn !! sign to be applied to acceleration end subroutine symba_kick_encounter - module subroutine symba_netcdf_initialize_output(self, param) - use swiftest_classes, only : swiftest_parameters - implicit none - class(symba_netcdf_parameters), intent(inout) :: self !! Parameters used to identify a particular NetCDF dataset - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters - end subroutine symba_netcdf_initialize_output - - module subroutine symba_netcdf_open(self, param) - use swiftest_classes, only : swiftest_parameters - implicit none - class(symba_netcdf_parameters), intent(inout) :: self !! Parameters used to identify a particular NetCDF dataset - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters - end subroutine symba_netcdf_open - - module subroutine symba_netcdf_write_frame_cb(self, iu, param) - use swiftest_classes, only : swiftest_parameters, netcdf_parameters - implicit none - class(symba_cb), intent(in) :: self !! Symba central body object - class(netcdf_parameters), intent(inout) :: iu !! Parameters used to identify a particular NetCDF dataset - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters - end subroutine symba_netcdf_write_frame_cb - - module subroutine symba_netcdf_write_frame_pl(self, iu, param) - use swiftest_classes, only : swiftest_parameters, netcdf_parameters - implicit none - class(symba_pl), intent(in) :: self !! Symba massive body object - class(netcdf_parameters), intent(inout) :: iu !! Parameters used to identify a particular NetCDF dataset - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters - end subroutine symba_netcdf_write_frame_pl - - module subroutine symba_netcdf_write_frame_tp(self, iu, param) - use swiftest_classes, only : swiftest_parameters, netcdf_parameters - implicit none - class(symba_tp), intent(in) :: self !! SyMBA test particle object - class(netcdf_parameters), intent(inout) :: iu !! Parameters used to identify a particular NetCDF dataset - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters - end subroutine symba_netcdf_write_frame_tp - - module subroutine symba_setup_initialize_particle_info_system(self, param) - use swiftest_classes, only : swiftest_parameters - implicit none - class(symba_nbody_system), intent(inout) :: self !! SyMBA nbody system object - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters with SyMBA extensions - end subroutine symba_setup_initialize_particle_info_system - module subroutine symba_setup_initialize_system(self, param) use swiftest_classes, only : swiftest_parameters implicit none diff --git a/src/netcdf/netcdf.f90 b/src/netcdf/netcdf.f90 index 274b05e3e..d2f7a0f80 100644 --- a/src/netcdf/netcdf.f90 +++ b/src/netcdf/netcdf.f90 @@ -120,6 +120,17 @@ module subroutine netcdf_initialize_output(self, param) call check( nf90_def_var(self%ncid, GMESCAPE_VARNAME, self%out_type, self%time_dimid, self%GMescape_varid) ) end if + if (self%ltrack_origin) then + call check( nf90_def_var(self%ncid, ORIGIN_TYPE_VARNAME, NF90_CHAR, [self%str_dimid, self%id_dimid], self%origin_type_varid) ) + call check( nf90_def_var(self%ncid, ORIGIN_TIME_VARNAME, self%out_type, self%id_dimid, self%origin_time_varid) ) + call check( nf90_def_var(self%ncid, ORIGIN_XHX_VARNAME, self%out_type, self%id_dimid, self%origin_xhx_varid) ) + call check( nf90_def_var(self%ncid, ORIGIN_XHY_VARNAME, self%out_type, self%id_dimid, self%origin_xhy_varid) ) + call check( nf90_def_var(self%ncid, ORIGIN_XHZ_VARNAME, self%out_type, self%id_dimid, self%origin_xhz_varid) ) + call check( nf90_def_var(self%ncid, ORIGIN_VHX_VARNAME, self%out_type, self%id_dimid, self%origin_vhx_varid) ) + call check( nf90_def_var(self%ncid, ORIGIN_VHY_VARNAME, self%out_type, self%id_dimid, self%origin_vhy_varid) ) + call check( nf90_def_var(self%ncid, ORIGIN_VHZ_VARNAME, self%out_type, self%id_dimid, self%origin_vhz_varid) ) + end if + return end subroutine netcdf_initialize_output @@ -197,6 +208,17 @@ module subroutine netcdf_open(self, param) call check( nf90_inq_varid(self%ncid, GMESCAPE_VARNAME, self%GMescape_varid) ) end if + if (self%ltrack_origin) then + call check( nf90_inq_varid(self%ncid, ORIGIN_TYPE_VARNAME, self%origin_type_varid)) + call check( nf90_inq_varid(self%ncid, ORIGIN_TIME_VARNAME, self%origin_time_varid)) + call check( nf90_inq_varid(self%ncid, ORIGIN_XHX_VARNAME, self%origin_xhx_varid)) + call check( nf90_inq_varid(self%ncid, ORIGIN_XHY_VARNAME, self%origin_xhy_varid)) + call check( nf90_inq_varid(self%ncid, ORIGIN_XHZ_VARNAME, self%origin_xhz_varid)) + call check( nf90_inq_varid(self%ncid, ORIGIN_VHX_VARNAME, self%origin_vhx_varid)) + call check( nf90_inq_varid(self%ncid, ORIGIN_VHY_VARNAME, self%origin_vhy_varid)) + call check( nf90_inq_varid(self%ncid, ORIGIN_VHZ_VARNAME, self%origin_vhz_varid)) + end if + return end subroutine netcdf_open @@ -219,7 +241,7 @@ module subroutine netcdf_write_frame_base(self, iu, param) tslot = int(param%ioutput, kind=I4B) + 1 select type(self) - class is (swiftest_body) + class is (swiftest_body) associate(n => self%nbody) if (n == 0) return allocate(ind(n)) @@ -245,6 +267,7 @@ module subroutine netcdf_write_frame_base(self, iu, param) call check( nf90_put_var(iu%ncid, iu%vhy_varid, self%vh(2, j), start=[idslot, tslot]) ) call check( nf90_put_var(iu%ncid, iu%vhz_varid, self%vh(3, j), start=[idslot, tslot]) ) end if + if ((param%out_form == EL) .or. (param%out_form == XVEL)) then call check( nf90_put_var(iu%ncid, iu%a_varid, self%a(j), start=[idslot, tslot]) ) call check( nf90_put_var(iu%ncid, iu%e_varid, self%e(j), start=[idslot, tslot]) ) @@ -253,27 +276,42 @@ module subroutine netcdf_write_frame_base(self, iu, param) call check( nf90_put_var(iu%ncid, iu%omega_varid, self%omega(j), start=[idslot, tslot]) ) call check( nf90_put_var(iu%ncid, iu%capm_varid, self%capm(j), start=[idslot, tslot]) ) end if - select type(pl => self) + + if (iu%ltrack_origin) then + charstring = trim(adjustl(self%info(j)%origin_type)) + strlen = len(charstring) + call check( nf90_put_var(iu%ncid, iu%origin_type_varid, charstring, start=[1, idslot], count=[strlen, 1]) ) + call check( nf90_put_var(iu%ncid, iu%origin_time_varid, self%info(j)%origin_time, start=[idslot]) ) + call check( nf90_put_var(iu%ncid, iu%origin_xhx_varid, self%info(j)%origin_xh(1), start=[idslot]) ) + call check( nf90_put_var(iu%ncid, iu%origin_xhy_varid, self%info(j)%origin_xh(2), start=[idslot]) ) + call check( nf90_put_var(iu%ncid, iu%origin_xhz_varid, self%info(j)%origin_xh(3), start=[idslot]) ) + call check( nf90_put_var(iu%ncid, iu%origin_vhx_varid, self%info(j)%origin_vh(1), start=[idslot]) ) + call check( nf90_put_var(iu%ncid, iu%origin_vhy_varid, self%info(j)%origin_vh(2), start=[idslot]) ) + call check( nf90_put_var(iu%ncid, iu%origin_vhz_varid, self%info(j)%origin_vh(3), start=[idslot]) ) + end if + + select type(self) class is (swiftest_pl) ! Additional output if the passed polymorphic object is a massive body - call check( nf90_put_var(iu%ncid, iu%Gmass_varid, pl%Gmass(j), start=[idslot, tslot]) ) + call check( nf90_put_var(iu%ncid, iu%Gmass_varid, self%Gmass(j), start=[idslot, tslot]) ) if (param%lrhill_present) then - call check( nf90_put_var(iu%ncid, iu%rhill_varid, pl%rhill(j), start=[idslot, tslot]) ) + call check( nf90_put_var(iu%ncid, iu%rhill_varid, self%rhill(j), start=[idslot, tslot]) ) end if if (param%lclose) then - call check( nf90_put_var(iu%ncid, iu%radius_varid, pl%radius(j), start=[idslot, tslot]) ) + call check( nf90_put_var(iu%ncid, iu%radius_varid, self%radius(j), start=[idslot, tslot]) ) end if if (param%lrotation) then - call check( nf90_put_var(iu%ncid, iu%Ip1_varid, pl%Ip(1, j), start=[idslot, tslot]) ) - call check( nf90_put_var(iu%ncid, iu%Ip2_varid, pl%Ip(2, j), start=[idslot, tslot]) ) - call check( nf90_put_var(iu%ncid, iu%Ip3_varid, pl%Ip(3, j), start=[idslot, tslot]) ) - call check( nf90_put_var(iu%ncid, iu%rotx_varid, pl%rot(1, j), start=[idslot, tslot]) ) - call check( nf90_put_var(iu%ncid, iu%roty_varid, pl%rot(2, j), start=[idslot, tslot]) ) - call check( nf90_put_var(iu%ncid, iu%rotz_varid, pl%rot(3, j), start=[idslot, tslot]) ) + call check( nf90_put_var(iu%ncid, iu%Ip1_varid, self%Ip(1, j), start=[idslot, tslot]) ) + call check( nf90_put_var(iu%ncid, iu%Ip2_varid, self%Ip(2, j), start=[idslot, tslot]) ) + call check( nf90_put_var(iu%ncid, iu%Ip3_varid, self%Ip(3, j), start=[idslot, tslot]) ) + call check( nf90_put_var(iu%ncid, iu%rotx_varid, self%rot(1, j), start=[idslot, tslot]) ) + call check( nf90_put_var(iu%ncid, iu%roty_varid, self%rot(2, j), start=[idslot, tslot]) ) + call check( nf90_put_var(iu%ncid, iu%rotz_varid, self%rot(3, j), start=[idslot, tslot]) ) end if if (param%ltides) then - call check( nf90_put_var(iu%ncid, iu%k2_varid, pl%k2(j), start=[idslot, tslot]) ) - call check( nf90_put_var(iu%ncid, iu%Q_varid, pl%Q(j), start=[idslot, tslot]) ) + call check( nf90_put_var(iu%ncid, iu%k2_varid, self%k2(j), start=[idslot, tslot]) ) + call check( nf90_put_var(iu%ncid, iu%Q_varid, self%Q(j), start=[idslot, tslot]) ) end if + end select end do end associate @@ -303,11 +341,25 @@ module subroutine netcdf_write_frame_base(self, iu, param) call check( nf90_put_var(iu%ncid, iu%k2_varid, self%k2, start=[idslot, tslot]) ) call check( nf90_put_var(iu%ncid, iu%Q_varid, self%Q, start=[idslot, tslot]) ) end if + + if (iu%ltrack_origin) then + charstring = trim(adjustl(self%info%origin_type)) + strlen = len(charstring) + call check( nf90_put_var(iu%ncid, iu%origin_type_varid, charstring, start=[1, idslot], count=[strlen, 1]) ) + call check( nf90_put_var(iu%ncid, iu%origin_time_varid, self%info%origin_time, start=[idslot]) ) + call check( nf90_put_var(iu%ncid, iu%origin_xhx_varid, self%info%origin_xh(1), start=[idslot]) ) + call check( nf90_put_var(iu%ncid, iu%origin_xhy_varid, self%info%origin_xh(2), start=[idslot]) ) + call check( nf90_put_var(iu%ncid, iu%origin_xhz_varid, self%info%origin_xh(3), start=[idslot]) ) + call check( nf90_put_var(iu%ncid, iu%origin_vhx_varid, self%info%origin_vh(1), start=[idslot]) ) + call check( nf90_put_var(iu%ncid, iu%origin_vhy_varid, self%info%origin_vh(2), start=[idslot]) ) + call check( nf90_put_var(iu%ncid, iu%origin_vhz_varid, self%info%origin_vh(3), start=[idslot]) ) + end if end select return end subroutine netcdf_write_frame_base + module subroutine netcdf_write_hdr_system(self, iu, param) !! author: David A. Minton !! diff --git a/src/setup/setup.f90 b/src/setup/setup.f90 index 02a44f091..c87aea559 100644 --- a/src/setup/setup.f90 +++ b/src/setup/setup.f90 @@ -67,15 +67,6 @@ module subroutine setup_construct_system(system, param) call util_exit(FAILURE) end select - select type(system) - class is (symba_nbody_system) - if (.not.allocated(system%cb%info)) allocate(symba_particle_info :: system%cb%info) - if (.not.allocated(param%nciu)) allocate(symba_netcdf_parameters :: param%nciu) - class default - if (.not.allocated(system%cb%info)) allocate(swiftest_particle_info :: system%cb%info) - if (.not.allocated(param%nciu)) allocate(netcdf_parameters :: param%nciu) - end select - return end subroutine setup_construct_system @@ -144,9 +135,45 @@ module subroutine setup_initialize_particle_info_system(self, param) implicit none ! Arguments class(swiftest_nbody_system), intent(inout) :: self !! Swiftest nbody system object - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters + class(swiftest_parameters), intent(inout) :: param !! Current run configuration parameters + ! Internals + logical :: ltrack_origin + integer(I4B) :: i associate(cb => self%cb, pl => self%pl, npl => self%pl%nbody, tp => self%tp, ntp => self%tp%nbody) + select type(param) + class is (symba_parameters) + param%nciu%ltrack_origin = param%lfragmentation + class default + param%nciu%ltrack_origin = .false. + end select + + select type(param) + class is (symba_parameters) + ltrack_origin = param%lfragmentation + class default + ltrack_origin = .false. + end select + + if (ltrack_origin) then + cb%info%origin_type = "Initial conditions" + cb%info%origin_time = param%t0 + cb%info%origin_xh(:) = 0.0_DP + cb%info%origin_vh(:) = 0.0_DP + do i = 1, self%pl%nbody + pl%info(i)%origin_type = "Initial conditions" + pl%info(i)%origin_time = param%t0 + pl%info(i)%origin_xh(:) = self%pl%xh(:,i) + pl%info(i)%origin_vh(:) = self%pl%vh(:,i) + end do + do i = 1, self%tp%nbody + tp%info(i)%origin_type = "Initial conditions" + tp%info(i)%origin_time = param%t0 + tp%info(i)%origin_xh(:) = self%tp%xh(:,i) + tp%info(i)%origin_vh(:) = self%tp%vh(:,i) + end do + end if + cb%info%particle_type = CB_TYPE_NAME call cb%dump_particle_info(param) if (npl > 0) then @@ -157,6 +184,7 @@ module subroutine setup_initialize_particle_info_system(self, param) tp%info(1:ntp)%particle_type = TP_TYPE_NAME call tp%dump_particle_info(param) end if + end associate return diff --git a/src/symba/symba_collision.f90 b/src/symba/symba_collision.f90 index b4505741b..3f5ebadc3 100644 --- a/src/symba/symba_collision.f90 +++ b/src/symba/symba_collision.f90 @@ -839,120 +839,109 @@ subroutine symba_collision_mergeaddsub(system, param, family, id_frag, Ip_frag, class is (symba_pl) select type(pl_discards => system%pl_discards) class is (symba_merger) - select type(info => pl%info) - class is (symba_particle_info) - associate(pl_adds => system%pl_adds, cb => system%cb) - - ! Add the family bodies to the subtraction list - nfamily = size(family(:)) - nfrag = size(m_frag(:)) - lmask(:) = .false. - lmask(family(:)) = .true. - pl%status(family(:)) = MERGED - nstart = pl_discards%nbody + 1 - nend = pl_discards%nbody + nfamily - call pl_discards%append(pl, lmask) - pl%ldiscard(family(:)) = .true. - pl%lcollision(family(:)) = .true. - - ! Record how many bodies were subtracted in this event - pl_discards%ncomp(nstart:nend) = nfamily - - ! Setup new bodies - allocate(plnew, mold=pl) - call plnew%setup(nfrag, param) - ibiggest = family(maxloc(pl%Gmass(family(:)), dim=1)) + associate(info => pl%info, pl_adds => system%pl_adds, cb => system%cb) + ! Add the family bodies to the subtraction list + nfamily = size(family(:)) + nfrag = size(m_frag(:)) + lmask(:) = .false. + lmask(family(:)) = .true. + pl%status(family(:)) = MERGED + nstart = pl_discards%nbody + 1 + nend = pl_discards%nbody + nfamily + call pl_discards%append(pl, lmask) + pl%ldiscard(family(:)) = .true. + pl%lcollision(family(:)) = .true. - ! Copy over identification, information, and physical properties of the new bodies from the fragment list - plnew%id(1:nfrag) = id_frag(1:nfrag) - param%maxid = param%maxid + nfrag - plnew%xb(:, 1:nfrag) = xb_frag(:, 1:nfrag) - plnew%vb(:, 1:nfrag) = vb_frag(:, 1:nfrag) + ! Record how many bodies were subtracted in this event + pl_discards%ncomp(nstart:nend) = nfamily + + ! Setup new bodies + allocate(plnew, mold=pl) + call plnew%setup(nfrag, param) + ibiggest = family(maxloc(pl%Gmass(family(:)), dim=1)) + + ! Copy over identification, information, and physical properties of the new bodies from the fragment list + plnew%id(1:nfrag) = id_frag(1:nfrag) + param%maxid = param%maxid + nfrag + plnew%xb(:, 1:nfrag) = xb_frag(:, 1:nfrag) + plnew%vb(:, 1:nfrag) = vb_frag(:, 1:nfrag) + do i = 1, nfrag + plnew%xh(:,i) = xb_frag(:, i) - cb%xb(:) + plnew%vh(:,i) = vb_frag(:, i) - cb%vb(:) + end do + plnew%mass(1:nfrag) = m_frag(1:nfrag) + plnew%Gmass(1:nfrag) = param%GU * m_frag(1:nfrag) + plnew%radius(1:nfrag) = rad_frag(1:nfrag) + plnew%density(1:nfrag) = m_frag(1:nfrag) / rad_frag(1:nfrag) + + select case(status) + case(DISRUPTION) + plnew%info(1:nfrag)%origin_type = "Disruption" + plnew%status(1:nfrag) = NEW_PARTICLE + plnew%info(1:nfrag)%origin_time = param%t + do i = 1, nfrag + write(info(i)%name, FRAGFMT) id_frag(i) + plnew%info(i)%origin_xh(:) = plnew%xh(:,i) + plnew%info(i)%origin_vh(:) = plnew%vh(:,i) + end do + case(SUPERCATASTROPHIC) + plnew%info(1:nfrag)%origin_type = "Supercatastrophic" + plnew%status(1:nfrag) = NEW_PARTICLE + plnew%info(1:nfrag)%origin_time = param%t do i = 1, nfrag - plnew%xh(:,i) = xb_frag(:, i) - cb%xb(:) - plnew%vh(:,i) = vb_frag(:, i) - cb%vb(:) + write(info(i)%name, FRAGFMT) id_frag(i) + plnew%info(i)%origin_xh(:) = plnew%xh(:,i) + plnew%info(i)%origin_vh(:) = plnew%vh(:,i) + end do + case(HIT_AND_RUN_DISRUPT) + plnew%info(1)%name = pl%info(ibiggest)%name + plnew%info(1)%origin_type = pl%info(ibiggest)%origin_type + plnew%info(1)%origin_xh(:) = pl%info(ibiggest)%origin_xh(:) + plnew%info(1)%origin_vh(:) = pl%info(ibiggest)%origin_vh(:) + plnew%status(1) = OLD_PARTICLE + plnew%status(2:nfrag) = NEW_PARTICLE + plnew%info(2:nfrag)%origin_type = "Hit and run fragment" + plnew%info(2:nfrag)%origin_time = param%t + do i = 2, nfrag + write(info(i)%name, FRAGFMT) id_frag(i) + plnew%info(i)%origin_xh(:) = plnew%xh(:,i) + plnew%info(i)%origin_vh(:) = plnew%vh(:,i) end do - plnew%mass(1:nfrag) = m_frag(1:nfrag) - plnew%Gmass(1:nfrag) = param%GU * m_frag(1:nfrag) - plnew%radius(1:nfrag) = rad_frag(1:nfrag) - plnew%density(1:nfrag) = m_frag(1:nfrag) / rad_frag(1:nfrag) + case(MERGED) + plnew%info(1)%name = pl%info(ibiggest)%name + plnew%info(1)%origin_type = pl%info(ibiggest)%origin_type + plnew%info(1)%origin_xh(:) = pl%info(ibiggest)%origin_xh(:) + plnew%info(1)%origin_vh(:) = pl%info(ibiggest)%origin_vh(:) + plnew%status(1) = OLD_PARTICLE + end select - select type(info => plnew%info) - class is (symba_particle_info) - select case(status) - case(DISRUPTION) - info(1:nfrag)%origin_type = "Disruption" - plnew%status(1:nfrag) = NEW_PARTICLE - info(1:nfrag)%origin_time = param%t - do i = 1, nfrag - write(info(i)%name, FRAGFMT) id_frag(i) - info(i)%origin_xh(:) = plnew%xh(:,i) - info(i)%origin_vh(:) = plnew%vh(:,i) - end do - case(SUPERCATASTROPHIC) - info(1:nfrag)%origin_type = "Supercatastrophic" - plnew%status(1:nfrag) = NEW_PARTICLE - info(1:nfrag)%origin_time = param%t - do i = 1, nfrag - write(info(i)%name, FRAGFMT) id_frag(i) - info(i)%origin_xh(:) = plnew%xh(:,i) - info(i)%origin_vh(:) = plnew%vh(:,i) - end do - case(HIT_AND_RUN_DISRUPT) - select type(plinfo => pl%info) - class is (symba_particle_info) - info(1)%name = plinfo(ibiggest)%name - info(1)%origin_xh(:) = plinfo(ibiggest)%origin_xh(:) - info(1)%origin_vh(:) = plinfo(ibiggest)%origin_vh(:) - end select - plnew%status(1) = OLD_PARTICLE - plnew%status(2:nfrag) = NEW_PARTICLE - info(2:nfrag)%origin_type = "Hit and run fragment" - info(2:nfrag)%origin_time = param%t - do i = 2, nfrag - write(info(i)%name, FRAGFMT) id_frag(i) - info(i)%origin_xh(:) = plnew%xh(:,i) - info(i)%origin_vh(:) = plnew%vh(:,i) - end do - case(MERGED) - select type(plinfo => pl%info) - class is (symba_particle_info) - info(1)%name = plinfo(ibiggest)%name - info(1)%origin_xh(:) = plinfo(ibiggest)%origin_xh(:) - info(1)%origin_vh(:) = plinfo(ibiggest)%origin_vh(:) - end select - plnew%status(1) = OLD_PARTICLE - end select - end select - - if (param%lrotation) then - plnew%Ip(:, 1:nfrag) = Ip_frag(:, 1:nfrag) - plnew%rot(:, 1:nfrag) = rot_frag(:, 1:nfrag) - end if - - if (param%ltides) then - plnew%Q = pl%Q(ibiggest) - plnew%k2 = pl%k2(ibiggest) - plnew%tlag = pl%tlag(ibiggest) - end if + if (param%lrotation) then + plnew%Ip(:, 1:nfrag) = Ip_frag(:, 1:nfrag) + plnew%rot(:, 1:nfrag) = rot_frag(:, 1:nfrag) + end if + + if (param%ltides) then + plnew%Q = pl%Q(ibiggest) + plnew%k2 = pl%k2(ibiggest) + plnew%tlag = pl%tlag(ibiggest) + end if - call plnew%set_mu(cb) - !Copy over or set integration parameters for new bodies - plnew%lcollision(1:nfrag) = .false. - plnew%ldiscard(1:nfrag) = .false. - plnew%levelg(1:nfrag) = pl%levelg(ibiggest) - plnew%levelm(1:nfrag) = pl%levelm(ibiggest) - - ! Append the new merged body to the list and record how many we made - nstart = pl_adds%nbody + 1 - nend = pl_adds%nbody + plnew%nbody - call pl_adds%append(plnew, lsource_mask=[(.true., i=1, plnew%nbody)]) - pl_adds%ncomp(nstart:nend) = plnew%nbody - - call plnew%setup(0, param) - deallocate(plnew) - end associate - end select + call plnew%set_mu(cb) + !Copy over or set integration parameters for new bodies + plnew%lcollision(1:nfrag) = .false. + plnew%ldiscard(1:nfrag) = .false. + plnew%levelg(1:nfrag) = pl%levelg(ibiggest) + plnew%levelm(1:nfrag) = pl%levelm(ibiggest) + + ! Append the new merged body to the list and record how many we made + nstart = pl_adds%nbody + 1 + nend = pl_adds%nbody + plnew%nbody + call pl_adds%append(plnew, lsource_mask=[(.true., i=1, plnew%nbody)]) + pl_adds%ncomp(nstart:nend) = plnew%nbody + + call plnew%setup(0, param) + deallocate(plnew) + end associate end select end select diff --git a/src/symba/symba_io.f90 b/src/symba/symba_io.f90 index 281fe7cdb..8b64950e9 100644 --- a/src/symba/symba_io.f90 +++ b/src/symba/symba_io.f90 @@ -2,30 +2,6 @@ use swiftest contains - module subroutine symba_io_dump_particle_info(self, iu) - !! author: David A. Minton - !! - !! Reads in particle information object information from an open file unformatted file - implicit none - ! Arguments - class(symba_particle_info), intent(in) :: self !! Particle metadata information object - integer(I4B), intent(in) :: iu !! Open file unit number - ! Internals - character(STRMAX) :: errmsg - - write(iu, err = 667, iomsg = errmsg) self%origin_type - write(iu, err = 667, iomsg = errmsg) self%origin_time - write(iu, err = 667, iomsg = errmsg) self%origin_xh(:) - write(iu, err = 667, iomsg = errmsg) self%origin_vh(:) - - return - - 667 continue - write(*,*) "Error writing particle metadata information from file: " // trim(adjustl(errmsg)) - call util_exit(FAILURE) - end subroutine symba_io_dump_particle_info - - module subroutine symba_io_param_reader(self, unit, iotype, v_list, iostat, iomsg) !! author: The Purdue Swiftest Team - David A. Minton, Carlisle A. Wishard, Jennifer L.L. Pouplin, and Jacob R. Elliott !! @@ -203,31 +179,6 @@ module subroutine symba_io_param_writer(self, unit, iotype, v_list, iostat, ioms end subroutine symba_io_param_writer - module subroutine symba_io_read_in_particle_info(self, iu) - !! author: David A. Minton - !! - !! Reads in particle information object information from an open file unformatted file - implicit none - ! Arguments - class(symba_particle_info), intent(inout) :: self !! Particle metadata information object - integer(I4B), intent(in) :: iu !! Open file unit number - ! Internals - character(STRMAX) :: errmsg - - call io_read_in_particle_info(self, iu) - read(iu, err = 667, iomsg = errmsg) self%origin_type - read(iu, err = 667, iomsg = errmsg) self%origin_time - read(iu, err = 667, iomsg = errmsg) self%origin_xh(:) - read(iu, err = 667, iomsg = errmsg) self%origin_vh(:) - - return - - 667 continue - write(*,*) "Error reading particle metadata information from file: " // trim(adjustl(errmsg)) - call util_exit(FAILURE) - end subroutine symba_io_read_in_particle_info - - module subroutine symba_io_write_discard(self, param) implicit none class(symba_nbody_system), intent(inout) :: self !! SyMBA nbody system object diff --git a/src/symba/symba_netcdf.f90 b/src/symba/symba_netcdf.f90 deleted file mode 100644 index a3ef5cec0..000000000 --- a/src/symba/symba_netcdf.f90 +++ /dev/null @@ -1,206 +0,0 @@ -submodule (symba_classes) s_symba_netcdf - use swiftest - use netcdf -contains - - subroutine check(status) - !! author: Carlisle A. Wishard, Dana Singh, and David A. Minton - !! - !! Checks the status of all NetCDF operations to catch errors - implicit none - ! Arguments - integer, intent (in) :: status - - if(status /= nf90_noerr) then - write(*,*) trim(nf90_strerror(status)) - call util_exit(FAILURE) - end if - - return - end subroutine check - - module subroutine symba_netcdf_initialize_output(self, param) - !! author: Carlisle A. Wishard, Dana Singh, and David A. Minton - !! - !! Initialize a NetCDF file system and defines all variables. - implicit none - ! Arguments - class(symba_netcdf_parameters), intent(inout) :: self !! Parameters used to identify a particular NetCDF dataset - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters - - - call netcdf_initialize_output(self, param) - - ! Define the variables - call check( nf90_def_var(self%ncid, ORIGIN_TYPE_VARNAME, NF90_CHAR, [self%str_dimid, self%id_dimid], self%origin_type_varid) ) - call check( nf90_def_var(self%ncid, ORIGIN_TIME_VARNAME, self%out_type, self%id_dimid, self%origin_time_varid) ) - call check( nf90_def_var(self%ncid, ORIGIN_XHX_VARNAME, self%out_type, self%id_dimid, self%origin_xhx_varid) ) - call check( nf90_def_var(self%ncid, ORIGIN_XHY_VARNAME, self%out_type, self%id_dimid, self%origin_xhy_varid) ) - call check( nf90_def_var(self%ncid, ORIGIN_XHZ_VARNAME, self%out_type, self%id_dimid, self%origin_xhz_varid) ) - call check( nf90_def_var(self%ncid, ORIGIN_VHX_VARNAME, self%out_type, self%id_dimid, self%origin_vhx_varid) ) - call check( nf90_def_var(self%ncid, ORIGIN_VHY_VARNAME, self%out_type, self%id_dimid, self%origin_vhy_varid) ) - call check( nf90_def_var(self%ncid, ORIGIN_VHZ_VARNAME, self%out_type, self%id_dimid, self%origin_vhz_varid) ) - - return - - end subroutine symba_netcdf_initialize_output - - - module subroutine symba_netcdf_open(self, param) - !! author: Carlisle A. Wishard, Dana Singh, and David A. Minton - !! - !! Opens a NetCDF file and does the variable inquiries to activate variable ids - implicit none - ! Arguments - class(symba_netcdf_parameters), intent(inout) :: self !! Parameters used to identify a particular NetCDF dataset - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters - ! Internals - integer(I4B) :: old_mode - - call netcdf_open(self, param) - - call check( nf90_inq_varid(self%ncid, ORIGIN_TYPE_VARNAME, self%origin_type_varid)) - call check( nf90_inq_varid(self%ncid, ORIGIN_TIME_VARNAME, self%origin_time_varid)) - call check( nf90_inq_varid(self%ncid, ORIGIN_XHX_VARNAME, self%origin_xhx_varid)) - call check( nf90_inq_varid(self%ncid, ORIGIN_XHY_VARNAME, self%origin_xhy_varid)) - call check( nf90_inq_varid(self%ncid, ORIGIN_XHZ_VARNAME, self%origin_xhz_varid)) - call check( nf90_inq_varid(self%ncid, ORIGIN_VHX_VARNAME, self%origin_vhx_varid)) - call check( nf90_inq_varid(self%ncid, ORIGIN_VHY_VARNAME, self%origin_vhy_varid)) - call check( nf90_inq_varid(self%ncid, ORIGIN_VHZ_VARNAME, self%origin_vhz_varid)) - - return - end subroutine symba_netcdf_open - - module subroutine symba_netcdf_write_frame_cb(self, iu, param) - !! author: Carlisle A. Wishard, Dana Singh, and David A. Minton - !! - !! Write a frame of output of a SyMBA massive body data to the binary output file - !! Note: If outputting to orbital elements, but sure that the conversion is done prior to calling this method - implicit none - ! Arguments - class(symba_cb), intent(in) :: self !! SyMBA central body object - class(netcdf_parameters), intent(inout) :: iu !! Parameters used to identify a particular NetCDF dataset - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters - ! Internals - integer(I4B) :: strlen, idslot - character(len=:), allocatable :: charstring - - call netcdf_write_frame_base(self, iu, param) - select type(iu) - class is (symba_netcdf_parameters) - select type(info => self%info) - class is (symba_particle_info) - idslot = self%id + 1 - - charstring = trim(adjustl(info%origin_type)) - strlen = len(charstring) - call check( nf90_put_var(iu%ncid, iu%origin_type_varid, charstring, start=[1, idslot], count=[strlen, 1]) ) - call check( nf90_put_var(iu%ncid, iu%origin_time_varid, info%origin_time, start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_xhx_varid, info%origin_xh(1), start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_xhy_varid, info%origin_xh(2), start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_xhz_varid, info%origin_xh(3), start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_vhx_varid, info%origin_vh(1), start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_vhy_varid, info%origin_vh(2), start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_vhz_varid, info%origin_vh(3), start=[idslot]) ) - end select - end select - - return - end subroutine symba_netcdf_write_frame_cb - - - module subroutine symba_netcdf_write_frame_pl(self, iu, param) - !! author: Carlisle A. Wishard, Dana Singh, and David A. Minton - !! - !! Write a frame of output of a SyMBA massive body data to the binary output file - !! Note: If outputting to orbital elements, but sure that the conversion is done prior to calling this method - implicit none - ! Arguments - class(symba_pl), intent(in) :: self !! SyMBA massive body object - class(netcdf_parameters), intent(inout) :: iu !! Parameters used to identify a particular NetCDF dataset - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters - ! Internals - integer(I4B) :: i, j, strlen, idslot - integer(I4B), dimension(:), allocatable :: ind - character(len=:), allocatable :: charstring - - call netcdf_write_frame_base(self, iu, param) - select type(iu) - class is (symba_netcdf_parameters) - associate(npl => self%nbody) - if (npl == 0) return - allocate(ind(npl)) - call util_sort(self%id(1:npl), ind(1:npl)) - select type(info => self%info) - class is (symba_particle_info) - do i = 1, npl - j = ind(i) - idslot = self%id(j) + 1 - - charstring = trim(adjustl(info(j)%origin_type)) - strlen = len(charstring) - call check( nf90_put_var(iu%ncid, iu%origin_type_varid, charstring, start=[1, idslot], count=[strlen, 1]) ) - call check( nf90_put_var(iu%ncid, iu%origin_time_varid, info(j)%origin_time, start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_xhx_varid, info(j)%origin_xh(1), start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_xhy_varid, info(j)%origin_xh(2), start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_xhz_varid, info(j)%origin_xh(3), start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_vhx_varid, info(j)%origin_vh(1), start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_vhy_varid, info(j)%origin_vh(2), start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_vhz_varid, info(j)%origin_vh(3), start=[idslot]) ) - end do - end select - end associate - - end select - - return - end subroutine symba_netcdf_write_frame_pl - - - module subroutine symba_netcdf_write_frame_tp(self, iu, param) - !! author: Carlisle A. Wishard, Dana Singh, and David A. Minton - !! - !! Write a frame of output of a SyMBA massive body data to the binary output file - !! Note: If outputting to orbital elements, but sure that the conversion is done prior to calling this method - implicit none - ! Arguments - class(symba_tp), intent(in) :: self !! SyMBA test particle - class(netcdf_parameters), intent(inout) :: iu !! Parameters used to identify a particular NetCDF dataset - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters - ! Internals - integer(I4B) :: i, j, strlen, idslot - integer(I4B), dimension(:), allocatable :: ind - character(len=:), allocatable :: charstring - - call netcdf_write_frame_base(self, iu, param) - select type(iu) - class is (symba_netcdf_parameters) - associate(ntp => self%nbody) - if (ntp == 0) return - allocate(ind(ntp)) - call util_sort(self%id(1:ntp), ind(1:ntp)) - select type(info => self%info) - class is (symba_particle_info) - do i = 1, ntp - j = ind(i) - idslot = self%id(j) + 1 - - charstring = trim(adjustl(info(j)%origin_type)) - strlen = len(charstring) - call check( nf90_put_var(iu%ncid, iu%origin_type_varid, charstring, start=[1, idslot], count=[strlen, 1]) ) - call check( nf90_put_var(iu%ncid, iu%origin_time_varid, info(j)%origin_time, start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_xhx_varid, info(j)%origin_xh(1), start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_xhy_varid, info(j)%origin_xh(2), start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_xhz_varid, info(j)%origin_xh(3), start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_vhx_varid, info(j)%origin_vh(1), start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_vhy_varid, info(j)%origin_vh(2), start=[idslot]) ) - call check( nf90_put_var(iu%ncid, iu%origin_vhz_varid, info(j)%origin_vh(3), start=[idslot]) ) - end do - end select - end associate - - end select - - return - end subroutine symba_netcdf_write_frame_tp -end submodule s_symba_netcdf \ No newline at end of file diff --git a/src/symba/symba_setup.f90 b/src/symba/symba_setup.f90 index f7cc82024..9eff6bf6c 100644 --- a/src/symba/symba_setup.f90 +++ b/src/symba/symba_setup.f90 @@ -2,50 +2,6 @@ use swiftest contains - module subroutine symba_setup_initialize_particle_info_system(self, param) - !! author: David A. Minton - !! - !! Initializes a new particle information data structure with initial conditions recorded - implicit none - ! Argumets - class(symba_nbody_system), intent(inout) :: self !! SyMBA nbody system object - class(swiftest_parameters), intent(in) :: param !! Current run configuration parameters - ! Internals - integer(I4B) :: i - - select type(cbinfo => self%cb%info) - class is (symba_particle_info) - cbinfo%origin_type = "Initial conditions" - cbinfo%origin_time = param%t0 - cbinfo%origin_xh(:) = 0.0_DP - cbinfo%origin_vh(:) = 0.0_DP - end select - - select type(plinfo => self%pl%info) - class is (symba_particle_info) - do i = 1, self%pl%nbody - plinfo(i)%origin_type = "Initial conditions" - plinfo(i)%origin_time = param%t0 - plinfo(i)%origin_xh(:) = self%pl%xh(:,i) - plinfo(i)%origin_vh(:) = self%pl%vh(:,i) - end do - end select - - select type(tpinfo => self%tp%info) - class is (symba_particle_info) - do i = 1, self%tp%nbody - tpinfo(i)%origin_type = "Initial conditions" - tpinfo(i)%origin_time = param%t0 - tpinfo(i)%origin_xh(:) = self%tp%xh(:,i) - tpinfo(i)%origin_vh(:) = self%tp%vh(:,i) - end do - end select - call setup_initialize_particle_info_system(self, param) - - return - end subroutine symba_setup_initialize_particle_info_system - - module subroutine symba_setup_initialize_system(self, param) !! author: David A. Minton !! @@ -116,7 +72,6 @@ module subroutine symba_setup_pl(self, n, param) if (n <= 0) return - if (allocated(self%info)) deallocate(self%info) if (allocated(self%lcollision)) deallocate(self%lcollision) if (allocated(self%lencounter)) deallocate(self%lencounter) if (allocated(self%lmtiny)) deallocate(self%lmtiny) @@ -128,9 +83,7 @@ module subroutine symba_setup_pl(self, n, param) if (allocated(self%peri)) deallocate(self%peri) if (allocated(self%atp)) deallocate(self%atp) if (allocated(self%kin)) deallocate(self%kin) - if (allocated(self%info)) deallocate(self%info) - allocate(symba_particle_info :: self%info(n)) allocate(self%lcollision(n)) allocate(self%lencounter(n)) allocate(self%lmtiny(n)) diff --git a/src/util/util_append.f90 b/src/util/util_append.f90 index 98a98ac0a..9722e4059 100644 --- a/src/util/util_append.f90 +++ b/src/util/util_append.f90 @@ -122,41 +122,23 @@ module subroutine util_append_arr_info(arr, source, nold, nsrc, lsource_mask) !! Append a single array of particle information type onto another. If the destination array is not allocated, or is not big enough, this will allocate space for it. implicit none ! Arguments - class(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Destination array - class(swiftest_particle_info), dimension(:), allocatable, intent(in) :: source !! Array to append - integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively - logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to + type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Destination array + type(swiftest_particle_info), dimension(:), allocatable, intent(in) :: source !! Array to append + integer(I4B), intent(in) :: nold, nsrc !! Extend of the old array and the source array, respectively + logical, dimension(:), intent(in) :: lsource_mask !! Logical mask indicating which elements to append to ! Internals integer(I4B) :: nnew - class(swiftest_particle_info), dimension(:), allocatable :: arr_tmp, source_tmp if (.not. allocated(source)) return nnew = count(lsource_mask(1:nsrc)) + if (.not.allocated(arr)) then + allocate(arr(nold+nnew)) + else + call util_resize(arr, nold + nnew) + end if - select type(source) - class is (symba_particle_info) - allocate(symba_particle_info :: arr_tmp(nold+nnew)) - if (nold > 0) then - arr_tmp(1:nold) = arr(1:nold) - deallocate(arr) - end if - class is (swiftest_particle_info) - allocate(swiftest_particle_info :: arr_tmp(nold+nnew)) - if (nold > 0) then - arr_tmp(1:nold) = arr(1:nold) - deallocate(arr) - end if - end select - - select type(source) - class is (symba_particle_info) - arr_tmp(nold + 1:nold + nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) - class is (swiftest_particle_info) - arr_tmp(nold + 1:nold + nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) - end select - - call move_alloc(arr_tmp, arr) + arr(nold + 1:nold + nnew) = pack(source(1:nsrc), lsource_mask(1:nsrc)) return end subroutine util_append_arr_info diff --git a/src/util/util_fill.f90 b/src/util/util_fill.f90 index f3f6a3a95..7009e3688 100644 --- a/src/util/util_fill.f90 +++ b/src/util/util_fill.f90 @@ -90,37 +90,18 @@ module subroutine util_fill_arr_info(keeps, inserts, lfill_list) !! This is the inverse of a spill operation implicit none ! Arguments - class(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - class(swiftest_particle_info), dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep + type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep + type(swiftest_particle_info), dimension(:), allocatable, intent(in) :: inserts !! Array of values to insert into keep logical, dimension(:), intent(in) :: lfill_list !! Logical array of bodies to merge into the keeps ! Internals - class(swiftest_particle_info), dimension(:), allocatable :: ktmp, itmp - integer(I4B) :: nk, ni + type(swiftest_particle_info), dimension(:), allocatable :: ktmp, itmp - if (.not.allocated(keeps) .or. .not.allocated(inserts)) return - - nk = size(keeps) - ni = size(inserts) - - select type(keeps) - class is (symba_particle_info) - allocate(symba_particle_info :: ktmp(nk)) - class is (swiftest_particle_info) - allocate(swiftest_particle_info :: ktmp(nk)) - end select - select type(inserts) - class is (symba_particle_info) - allocate(symba_particle_info :: itmp(ni)) - class is (swiftest_particle_info) - allocate(swiftest_particle_info :: itmp(ni)) - end select + if (.not.allocated(keeps) .or. .not.allocated(inserts)) return - ktmp(:) = unpack(ktmp(:), .not.lfill_list(:), ktmp(:)) - ktmp(:) = unpack(itmp(:), lfill_list(:), ktmp(:)) + keeps(:) = unpack(keeps(:), .not.lfill_list(:), keeps(:)) + keeps(:) = unpack(inserts(:), lfill_list(:), keeps(:)) - keeps(:) = ktmp(:) - return end subroutine util_fill_arr_info diff --git a/src/util/util_resize.f90 b/src/util/util_resize.f90 index 7ba7e19df..d43d0b879 100644 --- a/src/util/util_resize.f90 +++ b/src/util/util_resize.f90 @@ -138,28 +138,21 @@ module subroutine util_resize_arr_I4B(arr, nnew) end subroutine util_resize_arr_I4B - module subroutine util_resize_arr_info(arr, nnew) !! author: David A. Minton !! !! Resizes an array component of type character string. Array will only be resized if has previously been allocated. Passing nnew = 0 will deallocate. implicit none ! Arguments - class(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Array to resize + type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Array to resize integer(I4B), intent(in) :: nnew !! New size ! Internals - class(swiftest_particle_info), dimension(:), allocatable :: tmp !! Temporary storage array in case the input array is already allocated + type(swiftest_particle_info), dimension(:), allocatable :: tmp !! Temporary storage array in case the input array is already allocated integer(I4B) :: nold !! Old size logical :: is_symba if (.not. allocated(arr) .or. nnew < 0) return - select type(arr) - class is (symba_particle_info) - is_symba = .true. - class default - is_symba = .false. - end select nold = size(arr) if (nnew == nold) return @@ -167,12 +160,8 @@ module subroutine util_resize_arr_info(arr, nnew) deallocate(arr) return end if - - if (is_symba) then - allocate(symba_particle_info :: tmp(nnew)) - else - allocate(swiftest_particle_info :: tmp(nnew)) - end if + + allocate(tmp(nnew)) if (nnew > nold) then tmp(1:nold) = arr(1:nold) else diff --git a/src/util/util_sort.f90 b/src/util/util_sort.f90 index 4364aa02a..218e49ef3 100644 --- a/src/util/util_sort.f90 +++ b/src/util/util_sort.f90 @@ -465,22 +465,18 @@ module subroutine util_sort_rearrange_arr_info(arr, ind, n) !! Rearrange a single array of particle information type in-place from an index list. implicit none ! Arguments - class(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Destination array + type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: arr !! Destination array integer(I4B), dimension(:), intent(in) :: ind !! Index to rearrange against integer(I4B), intent(in) :: n !! Number of elements in arr and ind to rearrange ! Internals - class(swiftest_particle_info), dimension(:), allocatable :: tmp !! Temporary copy of array used during rearrange operation + type(swiftest_particle_info), dimension(:), allocatable :: tmp !! Temporary copy of array used during rearrange operation integer(I4B) :: i + if (.not. allocated(arr) .or. n <= 0) return - select type(arr) - class is (symba_particle_info) - allocate(symba_particle_info :: tmp(n)) - class is (swiftest_particle_info) - allocate(swiftest_particle_info :: tmp(n)) - end select + allocate(tmp, mold=arr) tmp(1:n) = arr(ind(1:n)) - arr(1:n) = tmp(1:n) + call move_alloc(tmp, arr) return end subroutine util_sort_rearrange_arr_info diff --git a/src/util/util_spill.f90 b/src/util/util_spill.f90 index 66e5d22d8..84391f168 100644 --- a/src/util/util_spill.f90 +++ b/src/util/util_spill.f90 @@ -208,8 +208,8 @@ module subroutine util_spill_arr_info(keeps, discards, lspill_list, ldestructive !! This is the inverse of a spill operation implicit none ! Arguments - class(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep - class(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: discards !! Array of discards + type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: keeps !! Array of values to keep + type(swiftest_particle_info), dimension(:), allocatable, intent(inout) :: discards !! Array of discards logical, dimension(:), intent(in) :: lspill_list !! Logical array of bodies to spill into the discardss logical, intent(in) :: ldestructive !! Logical flag indicating whether or not this operation should alter the keeps array or not ! Internals