Mercurial > repos > public > sbplib_julia

annotate EquidistantGrid.jl @ 95:9d53ecca34f7 cell_based_test

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Switch to using cartesian indicies
author Jonatan Werpers <jonatan@werpers.com>
date Mon, 04 Feb 2019 17:51:36 +0100
parents 8d505e9bc715
children 631eb9b35d72
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
51
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
1 # EquidistantGrid is a grid with equidistant grid spacing per coordinat
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
2 # direction. The domain is defined through the two points P1 = x̄₁, P2 = x̄₂
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
3 # by the exterior product of the vectors obtained by projecting (x̄₂-x̄₁) onto
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
4 # the coordinate directions. E.g for a 2D grid with x̄₁=(-1,0) and x̄₂=(1,2)
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
5 # the domain is defined as (-1,1)x(0,2).
75
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
6
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
7 struct EquidistantGrid{Dim,T<:Real} <: AbstractGrid
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
8 numberOfPointsPerDim::NTuple{Dim, Int} # First coordinate direction stored first, then
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
9
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
10 limit_lower::NTuple{Dim, T}
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
11 limit_upper::NTuple{Dim, T}
51
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
12
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
13 # General constructor
75
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
14 function EquidistantGrid(nPointsPerDim::NTuple{Dim, Int}, limit_lower::NTuple{Dim, T}, limit_upper::NTuple{Dim, T}) where Dim where T
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
15 @assert all(nPointsPerDim.>0)
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
16 @assert all(limit_upper.-limit_lower .!= 0)
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
17 return new{Dim,T}(nPointsPerDim, limit_lower, limit_upper)
51
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
18 end
75
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
19
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
20 # # 1D constructor which can be called as EquidistantGrid(m, (xl,xr))
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
21 # function EquidistantGrid(nPointsPerDim::Integer, lims::NTuple{2,Real})
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
22 # return EquidistantGrid((nPointsPerDim,), ((lims[1],),(lims[2],)))
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
23 # end
51
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
24
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
25 end
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
26
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
27 # Returns the number of dimensions of an EquidistantGrid.
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
28 #
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
29 # @Input: grid - an EquidistantGrid
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
30 # @Return: numberOfPoints - The number of dimensions
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
31 function numberOfDimensions(grid::EquidistantGrid)
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
32 return length(grid.numberOfPointsPerDim)
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
33 end
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
34
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
35 # Computes the total number of points of an EquidistantGrid.
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
36 #
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
37 # @Input: grid - an EquidistantGrid
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
38 # @Return: numberOfPoints - The total number of points
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
39 function numberOfPoints(grid::EquidistantGrid)
75
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
40 return prod(grid.numberOfPointsPerDim)
51
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
41 end
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
42
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
43 # Computes the grid spacing of an EquidistantGrid, i.e the unsigned distance
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
44 # between two points for each coordinate direction.
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
45 #
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
46 # @Input: grid - an EquidistantGrid
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
47 # @Return: h̄ - Grid spacing for each coordinate direction stored in a tuple.
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
48 function spacings(grid::EquidistantGrid)
75
93c833019857 Make EquidistantGrid more concrete
Jonatan Werpers <jonatan@werpers.com>
parents: 74
diff changeset
49 return abs.(grid.limit_upper.-grid.limit_lower)./(grid.numberOfPointsPerDim.-1)
51
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
50 end
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
51
95
9d53ecca34f7 Switch to using cartesian indicies
Jonatan Werpers <jonatan@werpers.com>
parents: 85
diff changeset
52 function Base.eachindex(grid::EquidistantGrid)
9d53ecca34f7 Switch to using cartesian indicies
Jonatan Werpers <jonatan@werpers.com>
parents: 85
diff changeset
53 CartesianIndices(grid.numberOfPointsPerDim)
9d53ecca34f7 Switch to using cartesian indicies
Jonatan Werpers <jonatan@werpers.com>
parents: 85
diff changeset
54 end
9d53ecca34f7 Switch to using cartesian indicies
Jonatan Werpers <jonatan@werpers.com>
parents: 85
diff changeset
55
51
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
56 # Computes the points of an EquidistantGrid as a vector of tuples. The vector is ordered
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
57 # such that points in the first coordinate direction varies first, then the second
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
58 # and lastely the third (if applicable)
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
59 #
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
60 # @Input: grid - an EquidistantGrid
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
61 # @Return: points - the points of the grid.
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
62 function points(grid::EquidistantGrid)
95
9d53ecca34f7 Switch to using cartesian indicies
Jonatan Werpers <jonatan@werpers.com>
parents: 85
diff changeset
63 # TODO: Make this return an abstract array?
9d53ecca34f7 Switch to using cartesian indicies
Jonatan Werpers <jonatan@werpers.com>
parents: 85
diff changeset
64 physical_domain_size = (grid.limit_upper .- grid.limit_lower)
9d53ecca34f7 Switch to using cartesian indicies
Jonatan Werpers <jonatan@werpers.com>
parents: 85
diff changeset
65 indices = Tuple.(CartesianIndices(grid.numberOfPointsPerDim))
9d53ecca34f7 Switch to using cartesian indicies
Jonatan Werpers <jonatan@werpers.com>
parents: 85
diff changeset
66 return broadcast(I -> grid.limit_lower .+ physical_domain_size.*(I.-1), indices)
51
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
67 end
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
68
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
69 function pointsalongdim(grid::EquidistantGrid, dim::Integer)
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
70 @assert dim<=numberOfDimensions(grid)
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
71 @assert dim>0
77
2be36b38389d Use limit_lower and limit_upper in pointsalongdim
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents: 75
diff changeset
72 points = range(grid.limit_lower[dim],stop=grid.limit_lower[dim],length=grid.numberOfPointsPerDim[dim])
51
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
73 end
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
74
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
75 using PyPlot, PyCall
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
76
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
77 function plotgridfunction(grid::EquidistantGrid, gridfunction)
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
78 if numberOfDimensions(grid) == 1
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
79 plot(pointsalongdim(grid,1), gridfunction, linewidth=2.0)
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
80 elseif numberOfDimensions(grid) == 2
58
f3ea15c8ebb6 Clean up plotgridfunction now that the matplotlib backend works. The issue was not in the julia code, but rather the installation of matplotlib.
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents: 53
diff changeset
81 mx = grid.numberOfPointsPerDim[1]
f3ea15c8ebb6 Clean up plotgridfunction now that the matplotlib backend works. The issue was not in the julia code, but rather the installation of matplotlib.
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents: 53
diff changeset
82 my = grid.numberOfPointsPerDim[2]
f3ea15c8ebb6 Clean up plotgridfunction now that the matplotlib backend works. The issue was not in the julia code, but rather the installation of matplotlib.
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents: 53
diff changeset
83 X = repeat(pointsalongdim(grid,1),1,my)
f3ea15c8ebb6 Clean up plotgridfunction now that the matplotlib backend works. The issue was not in the julia code, but rather the installation of matplotlib.
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents: 53
diff changeset
84 Y = permutedims(repeat(pointsalongdim(grid,2),1,mx))
f3ea15c8ebb6 Clean up plotgridfunction now that the matplotlib backend works. The issue was not in the julia code, but rather the installation of matplotlib.
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents: 53
diff changeset
85 plot_surface(X,Y,reshape(gridfunction,mx,my));
51
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
86 else
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
87 error(string("Plot not implemented for dimension ", string(numberOfDimensions(grid))))
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
88 end
614b56a017b9 Split grid.jl into AbstractGrid.jl and EquidistantGrid.jl
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
89 end