Calculate a variable sponge radius based on distance to the boundary
node's furthest neighbour.
(Adapted from Phil Hall's 'produce_netcdf_input_data.py')
spongeRadius = calc_sponge_radius(Mobj,Nlist)
DESCRIPTION
Calculates the sponge radius for each node on the open boundary, based
on the minimum of either the distance to the node's furthest
neighbour, or 100 km.
INPUT
Mobj = Matlab mesh object
Nlist = List of nodes
OUTPUT
spongeRadius = List of variable sponge radii
EXAMPLE USAGE
spongeRadius = calc_sponge_radius(Mobj,Nlist)
Author(s)
Karen Thurston (National Oceanography Centre, Liverpool)
==========================================================================0001 function [spongeRadius] = calc_sponge_radius(Mobj,Nlist) 0002 0003 % Calculate a variable sponge radius based on distance to the boundary 0004 % node's furthest neighbour. 0005 % (Adapted from Phil Hall's 'produce_netcdf_input_data.py') 0006 % 0007 % spongeRadius = calc_sponge_radius(Mobj,Nlist) 0008 % 0009 % DESCRIPTION 0010 % Calculates the sponge radius for each node on the open boundary, based 0011 % on the minimum of either the distance to the node's furthest 0012 % neighbour, or 100 km. 0013 % 0014 % INPUT 0015 % Mobj = Matlab mesh object 0016 % Nlist = List of nodes 0017 % 0018 % OUTPUT 0019 % spongeRadius = List of variable sponge radii 0020 % 0021 % EXAMPLE USAGE 0022 % spongeRadius = calc_sponge_radius(Mobj,Nlist) 0023 % 0024 % Author(s) 0025 % Karen Thurston (National Oceanography Centre, Liverpool) 0026 % 0027 %========================================================================== 0028 subname = 'calc_sponge_radius'; 0029 global ftbverbose 0030 if(ftbverbose) 0031 fprintf('\n') 0032 fprintf(['begin : ' subname '\n']) 0033 end 0034 0035 %-------------------------------------------------------------------------- 0036 % Get a unique list and make sure they are in the range of node numbers 0037 %-------------------------------------------------------------------------- 0038 Nlist = unique(Nlist); 0039 0040 spongeRadius = 100000+zeros(size(Nlist)); 0041 0042 % For each node on the open boundary 0043 for i =1:length(Nlist) 0044 % Find the neighbouring nodes 0045 [r,c]=find(Mobj.tri==Nlist(i)); 0046 neighbours = unique(Mobj.tri(r,:)); 0047 0048 % Remove the node of interest from the neighbours list 0049 n = find(neighbours~=Nlist(i)); 0050 neighbours = neighbours(n); 0051 0052 % Calculate the arc length (in degrees) between the node and its 0053 % neighbours 0054 arclen = distance(Mobj.lat(Nlist(i)),Mobj.lon(Nlist(i)),... 0055 Mobj.lat(neighbours),Mobj.lon(neighbours)); 0056 0057 % Convert from degrees to whole metres 0058 arclen = ceil(1000*deg2km(arclen)); 0059 0060 % If the smallest distance is less than 100km, keep it 0061 if min(arclen)<spongeRadius(i) 0062 spongeRadius(i)=min(arclen); 0063 end 0064 end 0065 0066 if(ftbverbose) 0067 fprintf(['end : ' subname '\n']) 0068 end 0069