changeset 1336:0666629aa183 feature/D2_boundary_opt

Add methods for creating grids with different grid point distributions for each coordinate direction, and also supports constructing periodic grids
author Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date Fri, 13 May 2022 13:26:16 +0200
parents 8d9fc7981796
children bf2554f1825d
files +grid/generalCurvilinear.m +multiblock/DefCurvilinear.m +util/get_periodic_grid.m
diffstat 3 files changed, 71 insertions(+), 8 deletions(-) [+]
line wrap: on
line diff
diff -r 8d9fc7981796 -r 0666629aa183 +grid/generalCurvilinear.m
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/+grid/generalCurvilinear.m	Fri May 13 13:26:16 2022 +0200
@@ -0,0 +1,54 @@
+% Creates a curvilinear grid of dimension length(m).
+% over the logical domain xi_lim, eta_lim, using a provided
+% coordinate mapping. The grid point distribution is
+% specified is 
+% Examples:
+%   g = grid.generalCurvilinear(mapping, [mx, my], xlim, ylim, {'equidist'}, {'boundaryopt', order, 'accurate'})
+function g = generalCurvilinear(mapping, m, varargin)
+    n = length(m);
+
+    % Check that parameters matches dimensions
+    matchingParams = false;
+    if length(varargin) == 2*n
+            matchingParams = iscell([varargin{1:2*n}]);
+    end
+    assert(matchingParams,'grid:generalCurvilinear:NonMatchingParameters','The number of parameters per dimensions do not match.');
+    
+    X = [];
+    h = [];
+    inds_periodic = [];
+    for i = 1:n
+        lim = varargin{i};
+        opts = varargin{i+n};
+        gridtype = opts{1};
+        switch gridtype
+        case 'equidist'
+            gridgenerator = @()util.get_grid(lim{1},lim{2},m(i));
+        case 'boundaryopt'
+            order = opts{2};
+            stencil_type = opts{3};
+            switch stencil_type
+            case {'Accurate','accurate','A','acc'}
+                gridgenerator = @()sbp.grid.accurateBoundaryOptimizedGrid(lim,m(i),order);
+            case {'Minimal','minimal','M','min'}
+                gridgenerator = @()sbp.grid.minimalBoundaryOptimizedGrid(lim,m(i),order);
+            end
+        case 'periodic'
+            gridgenerator = @()util.get_periodic_grid(lim{1},lim{2},m(i));
+            inds_periodic = [inds_periodic, i];
+        otherwise
+            error("grid type %s not supported. Must be one of 'equidist', 'boundaryopt', 'periodic'",gridtype);
+        end
+        try
+            [X{i},h(i)] = gridgenerator();
+        catch exception % Propagate any errors in the grid generation functions.
+            msgText = getReport(exception);
+            error('grid:boundaryOptimizedCurvilinear:InvalidParameter',msgText)
+        end
+    end
+    g = grid.Curvilinear(mapping, X{:});
+    g.logic.h = h;
+    for i = inds_periodic
+        g.logic.lim{i}{2} = g.logic.lim{i}{2}+h(i);
+    end
+end
\ No newline at end of file
diff -r 8d9fc7981796 -r 0666629aa183 +multiblock/DefCurvilinear.m
--- a/+multiblock/DefCurvilinear.m	Sat May 07 10:41:22 2022 +0200
+++ b/+multiblock/DefCurvilinear.m	Fri May 13 13:26:16 2022 +0200
@@ -52,10 +52,9 @@
             % If a scalar is passed, defer to getGridSizes implemented by subclass
             % TODO: This forces the interface of subclasses.
             % Should ms be included in varargin? Figure out bow to do it properly
-            if length(ms) == 1
+            if ~iscell(ms) && length(ms) == 1
                 ms = obj.getGridSizes(ms);
             end
-
             if isempty(varargin) || strcmp(varargin{1},'equidist')
                 gridgenerator = @(blockMap,m) grid.equidistantCurvilinear(blockMap, m);
             elseif strcmp(varargin{1},'boundaryopt')
@@ -63,13 +62,15 @@
                  stenciloption = varargin{3};
                  gridgenerator = @(blockMap,m) grid.boundaryOptimizedCurvilinear(blockMap,m,{0,1},{0,1},...
                      order,stenciloption);
-             else
-                 error('No grid type supplied!');
-             end
-             grids = cell(1, obj.nBlocks);
-             for i = 1:obj.nBlocks
+            elseif strcmp(varargin{1},'general')
+                gridgenerator = @(blockMap,m) grid.generalCurvilinear(blockMap,m,{0,1},{0,1},varargin{2:end});
+            else
+                error('No grid type supplied!');
+            end
+            grids = cell(1, obj.nBlocks);
+            for i = 1:obj.nBlocks
                 grids{i} = gridgenerator(obj.blockMaps{i}.S, ms{i});
-             end
+            end
 
             g = multiblock.Grid(grids, obj.connections, obj.boundaryGroups);
         end
diff -r 8d9fc7981796 -r 0666629aa183 +util/get_periodic_grid.m
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/+util/get_periodic_grid.m	Fri May 13 13:26:16 2022 +0200
@@ -0,0 +1,8 @@
+% Returns the 1D periodic grid where the last grid point at x_r
+% is omitted (it is the same as the the first grid point at x_l)
+function [x,h] = get_periodic_grid(x_l,x_r,m)
+    L = x_r-x_l;
+    h = L/m;
+    x = linspace(x_l,x_r,m+1)';
+    x = x(1:end-1);
+end
\ No newline at end of file