Granular.jl

Julia package for granular dynamics simulation
git clone git://src.adamsgaard.dk/Granular.jl
Log | Files | Refs | README | LICENSE

commit 12b44c7478f12efc5b9d145c56bdc32b66e1c43c
parent 3a3943c19320efec049abdfb35c442e0c58e6223
Author: Anders Damsgaard <andersd@riseup.net>
Date:   Thu, 22 Jun 2017 15:41:07 -0400

increase performance of cell coordinate functions by specificly passing coordinate arrays

Diffstat:
Msrc/datatypes.jl | 8++++----
Msrc/grid.jl | 30+++++++++++++++++-------------
Mtest/grid.jl | 4++--
3 files changed, 23 insertions(+), 19 deletions(-)

diff --git a/src/datatypes.jl b/src/datatypes.jl @@ -151,10 +151,10 @@ h-points. During read, the velocities are interpolated to the cell corners # Fields * `input_file::String`: path to input NetCDF file * `time::Array{Float64, 1}`: time in days -* `xq::Array{Float64, 1}`: nominal longitude of q-points [degrees_E] -* `yq::Array{Float64, 1}`: nominal latitude of q-points [degrees_N] -* `xh::Array{Float64, 1}`: nominal longitude of h-points [degrees_E] -* `yh::Array{Float64, 1}`: nominal latitude of h-points [degrees_N] +* `xq::Array{Float64, 2}`: nominal longitude of q-points [degrees_E] +* `yq::Array{Float64, 2}`: nominal latitude of q-points [degrees_N] +* `xh::Array{Float64, 2}`: nominal longitude of h-points [degrees_E] +* `yh::Array{Float64, 2}`: nominal latitude of h-points [degrees_N] * `zl::Array{Float64, 1}`: layer target potential density [kg m^-3] * `zi::Array{Float64, 1}`: interface target potential density [kg m^-3] * `u::Array{Float64, Int}`: zonal velocity (positive towards west) [m/s], diff --git a/src/grid.jl b/src/grid.jl @@ -67,7 +67,7 @@ function curl(grid::Any, k::Int, it::Int) - sw, se, ne, nw = getCellCornerCoordinates(grid, i, j) + sw, se, ne, nw = getCellCornerCoordinates(grid.xq, grid.yq, i, j) sw_se = norm(sw - se) se_ne = norm(se - ne) nw_ne = norm(nw - ne) @@ -209,7 +209,7 @@ This function is a wrapper for `getCellCornerCoordinates()` and function getNonDimensionalCellCoordinates(grid::Any, i::Int, j::Int, point::Array{float, 1}) - sw, se, ne, nw = getCellCornerCoordinates(grid, i, j) + sw, se, ne, nw = getCellCornerCoordinates(grid.xq, grid.yq, i, j) x_tilde, y_tilde = conformalQuadrilateralCoordinates(sw, se, ne, nw, point) return [x_tilde, y_tilde] end @@ -224,7 +224,7 @@ more robust. This function returns `true` or `false`. function isPointInCell(grid::Any, i::Int, j::Int, point::Array{float, 1}; method::String="Conformal") - sw, se, ne, nw = getCellCornerCoordinates(grid, i, j) + sw, se, ne, nw = getCellCornerCoordinates(grid.xq, grid.yq, i, j) if method == "Area" if areaOfQuadrilateral(sw, se, ne, nw) ≈ @@ -252,21 +252,23 @@ end export getCellCornerCoordinates """ - getCellCornerCoordinates(grid, i, j) + getCellCornerCoordinates(xq, yq, i, j) Returns grid corner coordinates in the following order (south-west corner, south-east corner, north-east corner, north-west corner). # Arguments -* `grid::Any`: grid object (Ocean or Atmosphere) containing grid. +* `xq::Array{Float64, 2}`: nominal longitude of q-points [degrees_E] +* `yq::Array{Float64, 2}`: nominal latitude of q-points [degrees_N] * `i::Int`: x-index of cell. * `j::Int`: y-index of cell. """ -function getCellCornerCoordinates(grid::Any, i::Int, j::Int) - sw = [grid.xq[ i, j], grid.yq[ i, j]] - se = [grid.xq[i+1, j], grid.yq[i+1, j]] - ne = [grid.xq[i+1, j+1], grid.yq[i+1, j+1]] - nw = [grid.xq[ i, j+1], grid.yq[ i, j+1]] +function getCellCornerCoordinates(xq::Array{Float64, 2}, yq::Array{Float64, 2}, + i::Int, j::Int) + sw = [xq[ i, j], yq[ i, j]] + se = [xq[i+1, j], yq[i+1, j]] + ne = [xq[i+1, j+1], yq[i+1, j+1]] + nw = [xq[ i, j+1], yq[ i, j+1]] return sw, se, ne, nw end @@ -277,12 +279,14 @@ export getCellCenterCoordinates Returns grid center coordinates (h-point). # Arguments -* `grid::Any`: grid object containing grid. +* `xh::Array{Float64, 2}`: nominal longitude of h-points [degrees_E] +* `yh::Array{Float64, 2}`: nominal latitude of h-points [degrees_N] * `i::Int`: x-index of cell. * `j::Int`: y-index of cell. """ -function getCellCenterCoordinates(grid::Any, i::Int, j::Int) - return [grid.xh[i, j], grid.yh[i, j]] +function getCellCenterCoordinates(xh::Array{Float64, 2}, yh::Array{Float64, 2}, + i::Int, j::Int) + return [xh[i, j], yh[i, j]] end export areaOfTriangle diff --git a/test/grid.jl b/test/grid.jl @@ -9,12 +9,12 @@ ocean = SeaIce.readOceanNetCDF("Baltic/00010101.ocean_month.nc", "Baltic/ocean_hgrid.nc") info("Testing coordinate retrieval functions") -sw, se, ne, nw = SeaIce.getCellCornerCoordinates(ocean, 1, 1) +sw, se, ne, nw = SeaIce.getCellCornerCoordinates(ocean.xq, ocean.yq, 1, 1) @test sw ≈ [6., 53.] @test se ≈ [7., 53.] @test ne ≈ [7., 54.] @test nw ≈ [6., 54.] -@test SeaIce.getCellCenterCoordinates(ocean, 1, 1) ≈ [6.5, 53.5] +@test SeaIce.getCellCenterCoordinates(ocean.xh, ocean.yh, 1, 1) ≈ [6.5, 53.5] info("Testing area-determination methods") @test SeaIce.areaOfTriangle([0., 0.], [1., 0.], [0., 1.]) ≈ .5