More user-friendly options for Laplacian solvers

[johnomotani]

Using flags to set many of the options for Laplacian solvers is cryptic, and not consistent with the way options are handled everywhere else. Could we replace them with something nicer? For example:

// Global options
bool zero_dc = false; /// Zero the DC (constant in Z) component of the solution
bool start_new = false; /// Iterative method start from solution=0. Has no effect for direct solvers
bool kx_zero = false; /// Zero the kx=0, n = 0 component

// Boundary options
std::string inner_boundary_condition; /// possible values "dirichlet", "neumann", "laplacian" or "cylinder"
std::string inner_boundary_condition_dc; /// default to same as inner_boundary_condition, possible values "dirichlet", "neumann", "laplacian", "cylinder", "gradpar", "gradparinv"
std::string inner_set; /// "zero", "initial_guess" or "rhs" [no flag, INVERT_SET or INVERT_RHS respectively]
bool inner_one_point = false; /// Only use one boundary point
... same for outer_* [except no "cylinder" option, as this is only needed for the inner boundary that could be the axis of the cylinder]

• Don't need replacement for INVERT_BOTH_BOUNDARY_ONE, can be set with boundary options
• Not sure INVERT_4th_ORDER option is needed - the serial_band solver is always 4th order anyway, the PETSc solver made an attempt at a 4th order implementation, but that's experimental/untested and I suspect the boundary conditions may not be correct (also the PETSc solver ignores the flag and reads a bool fourth_order option already)
• 'kx_zero'/INVERT_KX_ZERO is only supported by serial_tri and serial_band, maybe we can deprecate/remove?
• INVERT_SYM can be removed as it is not implemented anywhere (made redundant by re-definition of boundary conditions in 552c109)
• I think '*_one_point'/INVERT_BNDRY_ONE is not a sensible thing to have - if it is set then changing MXG changes the boundary condition. I'd suggest deprecating/removing this option. I don't know what it was originally for?

See also #91.

[ZedThree]

I am absolutely for this.

What's the best type for the *_boundary_condition options? Could users in principle define their own boundary conditions to be applied? If so, string makes sense. If they're things that needed to be implemented in the library, then an enum class might be better.

If lots of these options can be combined, bit flags might make sense, as opposed to taking lots of bools for instance.

[johnomotani]

I don't think there's a sensible way for users to define their own boundary conditions (without some major refactor that's beyond my imagination now anyway), so enum class probably makes most sense.

I don't think any of the options that have gone into *_boundary_condition/*_boundary_condition_dc can be combined, and I'm sure the values in *_set can't be.

I don't think the *_one_point option makes sense - for example it's inconsistent with allowing MXG to have different values because if MXG changes and *_one_point == true then the boundary condition applied in the Laplacian solver changes. So I'd be in favour of removing that option [updated list of suggestions in the post at the top to add this].

If we remove kx_zero and *_one_point there're only two bools left - that's not too many.

[johnomotani]

Just to make a note while I remember:
If we use an enum class (defined in the Laplacian namespace) like enum class BoundaryCondition {dirichlet, neumann, laplacian, cylinder, gradpar, gradparinv}; and include a function BoundaryCondition BoundaryConditionFromString(std::string s); then we can read the setting from a string option with something like

inner_boundary_condition = BoundaryConditionFromString(
                               options["inner_boundary_condition"]
                                   .withDefault("dirichlet")
                                   .doc("The inner boundary condition for the Laplacian solver."));

and this should throw an exception if the value passed is not supported, as long as we implement like e.g.

BOUT-dev/src/sys/bout_types.cxx

Lines 149 to 156 in 123e774

	YDirectionType YDirectionTypeFromString(const std::string& y_direction_string) {
	AUTO_TRACE();
	const static std::map<std::string, YDirectionType> stringToYDirectionType = {
	{"Standard", YDirectionType::Standard},
	{"Aligned", YDirectionType::Aligned}};
	return safeAt(stringToYDirectionType, y_direction_string);
	}

Another possibility would be to have a different enum class BoundaryCondition for each implementation, which would contain only the supported options, rather than having to check separately, but this probably is not so useful because the FFT-based solvers should all share the same set, which are implemented in tridagMatrix.