aboutsummaryrefslogtreecommitdiff
path: root/src/patch/patch_interp.cc
blob: 33358976358571aad87777925fb2bf399c4ea025 (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
// patch_frontier.cc -- interpolate gridfn data for another patch's ghost zone
// $Id$

//
// patch_frontier::patch_frontier
// patch_frontier::~patch_frontier
// patch_frontier::setup_interpolation_for_gridfn
// patch_frontier::interpolate
//

#include <stdio.h>
#include <assert.h>
#include <math.h>

#include "jt/stdc.h"
#include "jt/util.hh"
#include "jt/array.hh"
#include "jt/cpm_map.hh"
#include "jt/linear_map.hh"
#include "jt/interpolate.hh"

#include "fp.hh"
#include "coords.hh"
#include "grid.hh"
#include "fd_grid.hh"
#include "patch.hh"
#include "patch_edge.hh"
#include "ghost_zone.hh"
#include "patch_frontier.hh"

//*****************************************************************************
//*****************************************************************************
// patch_frontier - interpolation to get data for another patch's ghost zone
//*****************************************************************************
//*****************************************************************************

//
// This function constructs an  patch_frontier  object.
//
// Note this requires that the adjacent-side  ghost_zone  objects already
// exist, since they're used in determining the ipar range over which the
// interpolation will be done.
//
patch_frontier::patch_frontier(const patch_edge& my_edge_in,
			       int min_iperp_in, int max_iperp_in,
			       int min_parindex_in, int max_parindex_in,
			       const jtutil::array1d<fp>& min_parindex_used_in,
			       const jtutil::array1d<fp>& max_parindex_used_in,
			       const jtutil::array2d<fp>& interp_par_in,
			       int ghosted_min_gfn_in, int ghosted_max_gfn_in,
			       int interp_handle_in,
			       int interp_par_table_handle_in)
	: my_edge_(my_edge_in),
	  my_patch_(my_edge_in.my_patch()),
	  min_iperp_(min_iperp_in), max_iperp_(max_iperp_in),
	  ghosted_min_gfn_(ghosted_min_gfn_in),
	  ghosted_max_gfn_(ghosted_max_gfn_in),
	  interp_handle_(interp_handle_in),
	  interp_par_table_handle_(interp_par_table_handle_in),
	  interp_coord_origin_(0.0),	// computed later
	  interp_coord_delta_(my_patch().delta_ang(my_edge().par_is_rho())),
	  interp_par_coords_(ghosted_min_gfn_, ghosted_max_gfn_,
			     min_iperp_, max_iperp_,
			     min_parindex_in, max_parindex_in)
{
}

//*****************************************************************************

//
// This function destroys an  patch_frontier  object
//
patch_frontier::~patch_frontier()
{
// delete the interpolator objects
	for (int gfn = my_patch().min_gfn() ;
	     gfn <= my_patch().max_gfn() ;
	     ++gfn)
	{
	for (int iperp = min_iperp() ; iperp <= max_iperp() ; ++iperp)
	{
	delete (*interpolator_)(gfn,iperp);
	}
	}

delete interpolator_;
delete gridfn_buffer_;
delete ripar_map_;
}

//*****************************************************************************

//
// This function sets up a specified gridfn's interpolation state, i.e.
// it copies the gridfn data to the contiguous-in-ipar buffer and sets
// up any interpolator state
//
void patch_frontier::setup_interpolation_for_gridfn(int gfn)
{
	for (int iperp = min_iperp() ; iperp <= max_iperp() ; ++iperp)
	{
	//
	// copy this row of gridfn data from the (probably non-contiguous)
	// gridfn array, into our (contiguous-in-ipar) scratch array
	//
		for (int ipar = min_ipar() ; ipar <= max_ipar() ; ++ipar)
		{
		int irho   = my_edge().  irho_of_iperp_ipar(iperp, ipar);
		int isigma = my_edge().isigma_of_iperp_ipar(iperp, ipar);
		(*gridfn_buffer_)(gfn,iperp,ipar)
			= my_patch().gridfn(gfn,irho,isigma);
		}

	//
	// Now we want to set up the interpolator for the just-copied row
	// of data.  The interpolator interface requires that we pass a
	// pointer to the [0] data, where the [] subscripting refers to
	// the interpolator's integer coordinate.  Since ipar=0 may not
	// even exist, we introduce  ripar := ipar - min_ipar() and use
	// this as the interpolator's integer coordinate, so that [0]
	// refers to the start of this ipar row in the gridfn buffer.
	//
	const fp* row_ptr = & (*gridfn_buffer_)(gfn,iperp,min_ipar());
	interpolator_t_* I = (*interpolator_)(gfn,iperp);
	I->setup_data_y(row_ptr);
	}
}

//******************************************************************************

//
// This function does the actual interpolation of data for a specified
// gridfn to a specified (iperp,par) point.  It calls  error_exit()
// if the point is out of range or any other interpolator error is found.
//
fp patch_frontier::interpolate(int gfn, int iperp, fp par) const
{
const int N_dims = 1;