Skip to content

Instantly share code, notes, and snippets.

View briochemc's full-sized avatar

Benoît Pasquier briochemc

View GitHub Profile
using AIBECS, ModelingToolkit, LinearAlgebra
# dispatch for sparse matrices
Base.:*(A::SparseMatrixCSC{Operation,S}, x::StridedVector{Operation}) where {S} =
(T = Operation; mul!(similar(x, T, A.m), A, x, true, false))
Base.:*(A::SparseMatrixCSC{Tx,S}, x::StridedVector{Operation}) where {Tx,S} =
(T = LinearAlgebra.promote_op(LinearAlgebra.matprod, Operation, Tx); mul!(similar(x, T, A.m), A, x, true, false))
Base.:*(A::SparseMatrixCSC{Operation,S}, x::StridedVector{Tx}) where {Tx,S} =
(T = LinearAlgebra.promote_op(LinearAlgebra.matprod, Operation, Tx); mul!(similar(x, T, A.m), A, x, true, false))
using Luxor, Colors
function main(fname)
Drawing(600, 600, fname)
origin()
circle(O, 290, :clip)
sethue("white")
paint()
colors = [Luxor.julia_blue, Luxor.julia_red, Luxor.julia_green, Luxor.julia_purple]
using AIBECS, Plots
grd, TOCIM1 = OCIM1.load()
T(p) = TOCIM1
function G(x,p)
@unpack τ, z₀ = p
return @. 1 - x / τ * (z < z₀)
end
z = depthvec(grd)
F, ∇ₓF = state_function_and_Jacobian(T, G)
struct IdealAgeParameters{U} <: AbstractParameters{U}
# Load World Ocean Atlas 2D fields for PO₄, Si(OH)₄, and NO₄
using WorldOceanAtlasTools
product_year = 2018
tracer = "p"
period = 0
resolution = "1°"
ds = WorldOceanAtlasTools.WOA_Dataset(product_year, tracer, period, resolution, verbose=false) ;
field = "an"
field3D, lat, lon, depth = WorldOceanAtlasTools.get_gridded_3D_field(ds, tracer, field) ;
map_2D = field3D[:,:,1] ;
using SparseArrays, LinearAlgebra, SuiteSparse
using Unitful
"""
elunit(x)
Returns the unit of the elements of `x`.
(Returns a `Unitful.FreeUnits{(),NoDims,nothing}` object `x`'s unit is heterogeneous.)
"""
elunit(::SparseMatrixCSC{U,Int64}) where {U<:Quantity} = unit(U)
ENV["MPLBACKEND"]="qt5agg"
"qt5agg"
using PyPlot, PyCall
clf()
ax = subplot(projection=ccrs.NearsidePerspective(central_latitude=10.0, central_longitude=170.0, satellite_height=100000000.0))
ax.add_feature(cfeature.NaturalEarthFeature("physical", "ocean", "50m", edgecolor="face", facecolor="#4063d8"))
@briochemc
briochemc / resample.jl
Last active May 28, 2019 07:50
For prafter resampling script in Julia
length(ARGS) == 1 && (sheet = "Sheet1")
# Use the installed packages
using ExcelFiles, DataFrames, Interpolations
# Replace by your file with just the data!
df = DataFrame(load(abspath(ARGS[1]), sheet))
println("Resampling $(ARGS[1]) - $sheet")
using Makie
using GeometryTypes
wet3d = falses(2, 2, 2)
wet3d[[1, 2, 3, 5, 6]] .= true
shoebox3d = permutedims(wet3d, [2, 1, 3])
mycolor(x) = x ? RGBAf0(0,0,1,0.4) : RGBAf0(0.5,0.5,0,1)
Ind = CartesianIndices(size(shoebox3d))
using Makie
using AbstractPlotting
using FileIO, Colors
# Load World Ocean Atlas 2D fields for PO₄, Si(OH)₄, and NO₄
using WorldOceanAtlasTools
lat, lon, PO₄_2D = WorldOceanAtlasTools.WOA13_surface_map("p", 0, "1°")
lat, lon, NO₄_2D = WorldOceanAtlasTools.WOA13_surface_map("n", 0, "1°")
lat, lon, SiOH₄_2D = WorldOceanAtlasTools.WOA13_surface_map("i", 0, "1°")
# Add the image layer to each plot
# 1. Load packages
using DualMatrixTools, DualNumbers, SparseArrays
# 2. Create a random dual-valued vector x
n = 10
a, b = rand(n), rand(n) # real and dual parts
x = a + ε * b # dual-valued x
# 3. Create a sparse random dual-valued matrix M