### A Pluto.jl notebook ### # v0.19.40 using Markdown using InteractiveUtils # This Pluto notebook uses @bind for interactivity. When running this notebook outside of Pluto, the following 'mock version' of @bind gives bound variables a default value (instead of an error). macro bind(def, element) quote local iv = try Base.loaded_modules[Base.PkgId(Base.UUID("6e696c72-6542-2067-7265-42206c756150"), "AbstractPlutoDingetjes")].Bonds.initial_value catch; b -> missing; end local el = $(esc(element)) global $(esc(def)) = Core.applicable(Base.get, el) ? Base.get(el) : iv(el) el end end # ╔═╡ 104ce9b0-3fd1-11ec-3eff-3b029552e3d9 begin using IndividualDisplacements, CairoMakie, PlutoUI using OceanStateEstimation, MITgcmTools, CSV, JLD2 include("ECCO_FlowFields.jl") "Done with Loading Packages" end # ╔═╡ c9e9faa8-f5f0-479c-bc85-877ff7114883 md"""# Global Climatology Simulate the displacement of a set of individuals (or particles) according to climatological monthly mean flow at selected depth level (or in 3D). The flow fields are from a global ocean state estimate ([ECCO v4 r2](https://eccov4.readthedocs.io/en/latest/) ; see also ). It is here repeated for `ny` years to simulate trajectories. For additional documentation see : [1](https://JuliaClimate.github.io/IndividualDisplacements.jl/dev/), [2](https://JuliaClimate.github.io/MeshArrays.jl/dev/), [3](https://docs.juliadiffeq.org/latest/solvers/ode_solve.html), [4](https://en.wikipedia.org/wiki/Displacement_(vector)) [![simulated particle movie (5m)](https://user-images.githubusercontent.com/20276764/84766999-b801ad80-af9f-11ea-922a-610ad8a257dc.png)](https://youtu.be/W5DNqJG9jt0) """ # ╔═╡ 171fa252-7a35-4d4a-a940-60de77327cf4 TableOfContents() # ╔═╡ 7fec71b4-849f-4369-bec2-26bfe2e00a97 begin bind_k = (@bind ktxt Select(["0","1","10","30","40"],default="1")) bind_ny = (@bind nytxt Select(["1/12","1","2"],default="1/12")) bind_np = (@bind nptxt Select(["10","100","500"],default="10")) md"""## 1. Simulation Parameters The following parameters are used: - k = $(bind_k) = vertical level (for 2D; or 0 for 3D) - ny = $(bind_ny) = similated period (1/12 = 1 month) - np = $(bind_np) = number of individuals (100 by default) """ end # ╔═╡ 94ca10ae-6a8a-4038-ace0-07d7d9026712 md"""## 2. Configure `FlowFields` #### `global_ocean_circulation` reads in the global Ocean model grid and flow fields. It returns the `FlowFields` data structure `𝑃`. Ancillary variables are stored in NamedTuple `𝐷`. !!! note The global grid is handled via `MeshArrays.jl` and monthly climatology fields accessed via `OceanStateEstimation.jl`. """ # ╔═╡ 218b9beb-68f2-4498-a96d-08e0719b4cff begin k=parse(Int,ktxt) #vertical level or 0 np=parse(Int,nptxt) #number of individuals if nytxt=="1/12" ny=1 #number of years nm=1 #number of months else ny=parse(Int,nytxt) #number of years nm=12 #number of months end OceanStateEstimation.get_ecco_velocity_if_needed() 𝑃,𝐷=ECCO_FlowFields.global_ocean_circulation(k=k) "Done with Setting Up FlowFields" end # ╔═╡ f1215951-2eb2-490b-875a-91c1205b8f63 md"""## 3. Trajectory Computation ### 3.1 Initialization - initialize individual positions (`init_from_file`) - initial integration from time 0 to 0.5 month """ # ╔═╡ f727992f-b72a-45bc-93f1-cc8daf89af0f begin df = ECCO_FlowFields.init_from_file(np) 𝐼=Individuals(𝑃,df.x,df.y,df.f,(𝐷=𝐷,)) 📌_reference=deepcopy(𝐼.📌) 𝐼 end # ╔═╡ a3e45927-5d53-42be-b7b7-489d6e7a6fe5 begin 𝑇=(0.0,𝐼.𝑃.𝑇[2]) ∫!(𝐼,𝑇) ✔1="Done with Initial Integration" end # ╔═╡ 6158a5e4-89e0-4496-ab4a-044d1e3e8cc0 md""" ### 3.2 Monthly Step Function Since the Ocean circulation varies from month to month, the `𝐼.𝐷.🔄` function is used to update flow fields from one month to the next. Time variable flow fields are easily handled by defining a `step!` function that wraps around `∫!`. Here it does three things - `𝐼.𝐷.🔄` resets the velocity inputs (`𝐼.𝐷.u0` etc) to bracket time `t_ϵ=𝐼.𝑃.𝑇[2]+eps(𝐼.𝑃.𝑇[2])` - `reset_📌!` randomly selects a fraction of the individuals and resets their positions. This can be useful to maintain homogeneous coverage of the domain by the fleet of individuals. - `∫!(𝐼)` solves for the individual trajectories over one month (`𝐼.𝑃.𝑇`) with updated velocity fields (`𝐼.𝑃.u0` etc). It also adds diagnostics to the `DataFrame` used to record variables along the trajectory (`𝐼.🔴`). """ # ╔═╡ a2375720-f599-43b9-a7fb-af17956309b6 function step!(𝐼::Individuals) t_ϵ=𝐼.𝑃.𝑇[2]+eps(𝐼.𝑃.𝑇[2]) 𝐼.𝐷.🔄(𝐼.𝑃,𝐼.𝐷,t_ϵ) ECCO_FlowFields.reset_📌!(𝐼,𝐼.𝐷.frac,📌_reference) ∫!(𝐼) end # ╔═╡ 7efadea7-4542-40cf-893a-40a75e9c52be md"""### 3.3 Monthly Simulation Loop !!! note `add_lonlat!` derives geographic locations (longitude and latitude) from local grid coordinates (x, y, etc). """ # ╔═╡ 1044c5aa-1a56-45b6-a4c6-63d24eea878d begin ✔1 [step!(𝐼) for y=1:ny, m=1:nm] add_lonlat!(𝐼.🔴,𝐼.𝐷.XC,𝐼.𝐷.YC,𝑃.update_location!) ✔2="Done with Main Computation" end # ╔═╡ 6e43a2af-bf01-4f42-a4ba-1874a8cf4885 let ✔2 using DataFrames, Statistics gdf = groupby(𝐼.🔴, :ID) sgdf= combine(gdf,nrow,:lat => mean) end # ╔═╡ fc16b761-8b1f-41de-b4fe-7fa9987d6167 𝐼.🔴 # ╔═╡ c5ba37e9-2a68-4448-a2cb-dea1fbf08f1e md"""## 4. Visualize Displacements""" # ╔═╡ b4841dc0-c257-45e0-8657-79121f2c9ce8 plot(𝐼) # ╔═╡ 15077957-64d5-46a5-8a87-a76ad619cf38 md"""## 5. Summary Statistics Here we briefly demontrate the use of [DataFrames.jl](https://juliadata.github.io/DataFrames.jl/latest/) to analyze the output (𝐼.🔴) of our simulation. """ # ╔═╡ de8dbb43-68bc-4fb2-b0c8-07100b8a97a0 md"""## Appendix : Plotting Function""" # ╔═╡ e1cdcac9-c3cc-4ce4-a477-452ca460a3d5 """ plot(𝐼::Individuals) Plot initial and final positions, superimposed on a globalmap of ocean depth log. """ function plot(𝐼::Individuals) 𝐵=𝐼.𝐷.ODL xlims=extrema(𝐵.lon) ylims=extrema(𝐵.lat) fig=Figure() ax=Axis(fig[1,1]) plt=contourf!(ax,𝐵.lon,𝐵.lat,permutedims(𝐵.fld),colorrange=𝐵.rng,colormap = :ice) 🔴_by_t = groupby(𝐼.🔴, :t) lo=deepcopy(🔴_by_t[1].lon); lo[findall(lo.