Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Feature/fv3core fortran api #395

Merged
merged 22 commits into from
Dec 9, 2022
Merged
Show file tree
Hide file tree
Changes from 11 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions fv3core/pace/fv3core/__init__.py
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
from ._config import DynamicalCoreConfig
from .initialization.dycore_state import DycoreState
from .initialization.geos_wrapper import GeosDycoreWrapper
from .stencils.fv_dynamics import DynamicalCore
from .stencils.fv_subgridz import DryConvectiveAdjustment

Expand Down
256 changes: 256 additions & 0 deletions fv3core/pace/fv3core/initialization/geos_wrapper.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,256 @@
from datetime import timedelta
from typing import Dict

import f90nml
import numpy as np

import pace.util
from pace import fv3core


class GeosDycoreWrapper:
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

[nit] I believe NASA always refers it as GEOS

"""
Provides an interface for the Geos model to access the Pace dycore.
Takes numpy arrays as inputs, returns a dictionary of numpy arrays as outputs
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So this only works for numpy and cpu backend? We should extend this to take cupy array as well?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ideally we can also make this work on GPU. The reason for the current input/output setup is to comport with Purnendu's code, which puts the Fortran data in numpy arrays and reads a dictionary back out.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should be gpu-functional now, I hope?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This should be safe. Probably you want to had some comments in there to explain the upload mechanism.

Fortran comes as CPU
Dycore state is init two a GPU storage thanks to the backend in QuantityFactory
safe_assign_arrayis capable of recognizing that source & destination are not in the same memory space and deals with the upload / download.

In case Fortran memory is on GPU (unlikely) and comes as a cupy array this will take a copy

This upload mechanism is slow and requires it's own optimization via staggered cudaMemCpyAsync. This is not present at the moment

"""

def __init__(self, namelist: f90nml.Namelist, comm: pace.util.Comm, backend: str):
self.namelist = namelist

self.dycore_config = fv3core.DynamicalCoreConfig.from_f90nml(self.namelist)

self.layout = self.dycore_config.layout
partitioner = pace.util.CubedSpherePartitioner(
pace.util.TilePartitioner(self.layout)
)
self.communicator = pace.util.CubedSphereCommunicator(comm, partitioner)

sizer = pace.util.SubtileGridSizer.from_namelist(
self.namelist, partitioner.tile, self.communicator.tile.rank
)
quantity_factory = pace.util.QuantityFactory.from_backend(
sizer=sizer, backend=backend
)

# set up the metric terms and grid data
metric_terms = pace.util.grid.MetricTerms(
quantity_factory=quantity_factory, communicator=self.communicator
)
grid_data = pace.util.grid.GridData.new_from_metric_terms(metric_terms)

stencil_config = pace.dsl.stencil.StencilConfig(
compilation_config=pace.dsl.stencil.CompilationConfig(
backend=backend, rebuild=False, validate_args=False
),
)

self._grid_indexing = pace.dsl.stencil.GridIndexing.from_sizer_and_communicator(
sizer=sizer, cube=self.communicator
)
stencil_factory = pace.dsl.StencilFactory(
config=stencil_config, grid_indexing=self._grid_indexing
)

self.dycore_state = fv3core.DycoreState.init_zeros(
quantity_factory=quantity_factory
)

self.dycore_state.bdt = float(namelist["dt_atmos"])
if "fv_core_nml" in namelist.keys():
self.dycore_state.bdt = (
float(namelist["dt_atmos"]) / namelist["fv_core_nml"]["k_split"]
)
elif "dycore_config" in namelist.keys():
self.dycore_state.bdt = (
float(namelist["dt_atmos"]) / namelist["dycore_config"]["k_split"]
)
else:
raise KeyError("Cannot find k_split in namelist")

damping_coefficients = pace.util.grid.DampingCoefficients.new_from_metric_terms(
metric_terms
)

self.dynamical_core = fv3core.DynamicalCore(
comm=self.communicator,
grid_data=grid_data,
stencil_factory=stencil_factory,
quantity_factory=quantity_factory,
damping_coefficients=damping_coefficients,
config=self.dycore_config,
timestep=timedelta(seconds=self.dycore_config.dt_atmos),
phis=self.dycore_state.phis,
state=self.dycore_state,
)

self.output_dictionary: Dict[str, np.ndarray] = {}

def __call__(
self,
u: np.ndarray,
v: np.ndarray,
w: np.ndarray,
delz: np.ndarray,
pt: np.ndarray,
delp: np.ndarray,
q: np.ndarray,
ps: np.ndarray,
pe: np.ndarray,
pk: np.ndarray,
peln: np.ndarray,
pkz: np.ndarray,
phis: np.ndarray,
q_con: np.ndarray,
omga: np.ndarray,
ua: np.ndarray,
va: np.ndarray,
uc: np.ndarray,
vc: np.ndarray,
mfxd: np.ndarray,
mfyd: np.ndarray,
cxd: np.ndarray,
cyd: np.ndarray,
diss_estd: np.ndarray,
) -> dict:
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is there a reason for not making this more strongly typed? Could it be Dict[str, np.ndarray]?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

small oversight, thought I'd done that


self._put_fortran_data_in_dycore(
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Grab the timer from pace.util and slap one around each of those function call. See fv_dynamics for usage. They are GPU robust already

u,
v,
w,
delz,
pt,
delp,
q,
ps,
pe,
pk,
peln,
pkz,
phis,
q_con,
omga,
ua,
va,
uc,
vc,
mfxd,
mfyd,
cxd,
cyd,
diss_estd,
)

self.dynamical_core.step_dynamics(
state=self.dycore_state,
)

self._prep_outputs_for_geos()

return self.output_dictionary

def _put_fortran_data_in_dycore(
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Just wanted to check - is the idea here to not use translate fvdynamics anymore?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yup, that's the idea

self,
u: np.ndarray,
v: np.ndarray,
w: np.ndarray,
delz: np.ndarray,
pt: np.ndarray,
delp: np.ndarray,
q: np.ndarray,
ps: np.ndarray,
pe: np.ndarray,
pk: np.ndarray,
peln: np.ndarray,
pkz: np.ndarray,
phis: np.ndarray,
q_con: np.ndarray,
omga: np.ndarray,
ua: np.ndarray,
va: np.ndarray,
uc: np.ndarray,
vc: np.ndarray,
mfxd: np.ndarray,
mfyd: np.ndarray,
cxd: np.ndarray,
cyd: np.ndarray,
diss_estd: np.ndarray,
):
isc = self._grid_indexing.isc
jsc = self._grid_indexing.jsc
iec = self._grid_indexing.iec + 1
jec = self._grid_indexing.jec + 1

# Assign compute domain:
self.dycore_state.u.view[:] = u[isc:iec, jsc : jec + 1, :-1]
self.dycore_state.v.view[:] = v[isc : iec + 1, jsc:jec, :-1]
self.dycore_state.w.view[:] = w[isc:iec, jsc:jec, :-1]
self.dycore_state.ua.view[:] = ua[isc:iec, jsc:jec, :-1]
self.dycore_state.va.view[:] = va[isc:iec, jsc:jec, :-1]
self.dycore_state.uc.view[:] = uc[isc : iec + 1, jsc:jec, :-1]
self.dycore_state.vc.view[:] = vc[isc:iec, jsc : jec + 1, :-1]

self.dycore_state.delz.view[:] = delz[isc:iec, jsc:jec, :-1]
self.dycore_state.pt.view[:] = pt[isc:iec, jsc:jec, :-1]
self.dycore_state.delp.view[:] = delp[isc:iec, jsc:jec, :-1]

self.dycore_state.mfxd.view[:] = mfxd[isc : iec + 1, jsc:jec, :-1]
self.dycore_state.mfyd.view[:] = mfyd[isc:iec, jsc : jec + 1, :-1]
self.dycore_state.cxd.view[:] = cxd[isc : iec + 1, jsc:jec, :-1]
self.dycore_state.cyd.view[:] = cyd[isc:iec, jsc : jec + 1, :-1]

self.dycore_state.ps.view[:] = ps[isc:iec, jsc:jec]
self.dycore_state.pe.view[:] = pe[isc:iec, jsc:jec, :]
self.dycore_state.pk.view[:] = pk[isc:iec, jsc:jec, :]
self.dycore_state.peln.view[:] = peln[isc:iec, jsc:jec, :]
self.dycore_state.pkz.view[:] = pkz[isc:iec, jsc:jec, :-1]
self.dycore_state.phis.view[:] = phis[isc:iec, jsc:jec]
self.dycore_state.q_con.view[:] = q_con[isc:iec, jsc:jec, :-1]
self.dycore_state.omga.view[:] = omga[isc:iec, jsc:jec, :-1]
self.dycore_state.diss_estd.view[:] = diss_estd[isc:iec, jsc:jec, :-1]

# tracer quantities should be a 4d array in order:
# vapor, liquid, ice, rain, snow, graupel, cloud
self.dycore_state.qvapor.view[:] = q[isc:iec, jsc:jec, :-1, 0]
self.dycore_state.qliquid.view[:] = q[isc:iec, jsc:jec, :-1, 1]
self.dycore_state.qice.view[:] = q[isc:iec, jsc:jec, :-1, 2]
self.dycore_state.qrain.view[:] = q[isc:iec, jsc:jec, :-1, 3]
self.dycore_state.qsnow.view[:] = q[isc:iec, jsc:jec, :-1, 4]
self.dycore_state.qgraupel.view[:] = q[isc:iec, jsc:jec, :-1, 5]
if self.namelist["dycore_config"]["nwat"] > 6:
self.dycore_state.qcld.view[:] = q[isc:iec, jsc:jec, :-1, 6]

def _prep_outputs_for_geos(self):
self.output_dictionary["u"] = self.dycore_state.u.data[:]
self.output_dictionary["v"] = self.dycore_state.v.data[:]
self.output_dictionary["w"] = self.dycore_state.w.data[:]
self.output_dictionary["ua"] = self.dycore_state.ua.data[:]
self.output_dictionary["va"] = self.dycore_state.va.data[:]
self.output_dictionary["uc"] = self.dycore_state.uc.data[:]
self.output_dictionary["vc"] = self.dycore_state.vc.data[:]

self.output_dictionary["delc"] = self.dycore_state.delz.data[:]
self.output_dictionary["pt"] = self.dycore_state.pt.data[:]
self.output_dictionary["delp"] = self.dycore_state.delp.data[:]

self.output_dictionary["mfxd"] = self.dycore_state.mfxd.data[:]
self.output_dictionary["mfyd"] = self.dycore_state.mfyd.data[:]
self.output_dictionary["cxd"] = self.dycore_state.cxd.data[:]
self.output_dictionary["cyd"] = self.dycore_state.cyd.data[:]

self.output_dictionary["ps"] = self.dycore_state.ps.data[:]
self.output_dictionary["pe"] = self.dycore_state.pe.data[:]
self.output_dictionary["pk"] = self.dycore_state.pk.data[:]
self.output_dictionary["peln"] = self.dycore_state.peln.data[:]
self.output_dictionary["pkz"] = self.dycore_state.pkz.data[:]
self.output_dictionary["phis"] = self.dycore_state.phis.data[:]
self.output_dictionary["q_con"] = self.dycore_state.q_con.data[:]
self.output_dictionary["omga"] = self.dycore_state.omga.data[:]
self.output_dictionary["diss_estd"] = self.dycore_state.diss_estd.data[:]

self.output_dictionary["qvapor"] = self.dycore_state.qvapor.data[:]
self.output_dictionary["qliquid"] = self.dycore_state.qliquid.data[:]
self.output_dictionary["qice"] = self.dycore_state.qice.data[:]
self.output_dictionary["qrain"] = self.dycore_state.qrain.data[:]
self.output_dictionary["qsnow"] = self.dycore_state.qsnow.data[:]
self.output_dictionary["qgraupel"] = self.dycore_state.qgraupel.data[:]
self.output_dictionary["qcld"] = self.dycore_state.qcld.data[:]
Loading