Support complete elliptic integrals of 1st and 2nd type in Parser (AM…

…ReX-Codes#3225)

We want to use complete elliptic integrals from `cmath`, `comp_ellint_1`
and `comp_ellint_2` in the parser in order to make it easier to load in
magnetic fields from coils in WarpX. I added these to the parser, and
tested against the same functions in scipy - the results match.

Since clang does not support special functions in C++17
(https://en.cppreference.com/w/cpp/compiler_support/17#C.2B.2B17_library_features)
I implemented it for GCC compilers only, following the structure of how
`jn` is implemented in the parser.

Co-authored-by: Avigdor Veksler <aveksler@TAE7750-MLAP.tae.trialphaenergy.com>

Docs/sphinx_documentation/source/Basics.rst

@@ -399,8 +399,9 @@ expressions given in the form of string.  It supports ``+``, ``-``, ``*``,
 numbers can be computed with ``min`` and ``max``, respectively.  It supports
 the Heaviside step function, ``heaviside(x1,x2)`` that gives ``0``, ``x2``,
 ``1``, for ``x1 < 0``, ``x1 = 0`` and ``x1 > 0``, respectively.
-It also supports the Bessel function of the first kind of order ``n``
-``jn(n,x)``.
+It supports the Bessel function of the first kind of order ``n``
+``jn(n,x)``. Complete elliptic integrals of the first and second kind, ``comp_ellint_1`` and ``comp_ellint_2``,
+are supported only for gcc and CPUs.
 There is ``if(a,b,c)`` that gives ``b`` or ``c`` depending on the value of
 ``a``.  A number of comparison operators are supported, including ``<``,
 ``>``, ``==``, ``!=``, ``<=``, and ``>=``.  The Boolean results from

Src/Base/Parser/AMReX_Parser_Y.H

@@ -42,7 +42,9 @@ enum parser_f1_t {  // Built-in functions with one argument
     PARSER_POW_M1,
     PARSER_POW_P1,
     PARSER_POW_P2,
-    PARSER_POW_P3
+    PARSER_POW_P3,
+    PARSER_COMP_ELLINT_1,
+    PARSER_COMP_ELLINT_2
 };
 
 enum parser_f2_t {  // Built-in functions with two arguments
@@ -252,6 +254,32 @@ template <typename T>
 AMREX_GPU_HOST_DEVICE AMREX_NO_INLINE
 T parser_math_tanh (T a) { return std::tanh(a); }
 
+template <typename T>
+AMREX_GPU_HOST_DEVICE AMREX_NO_INLINE
+T parser_math_comp_ellint_1 (T a)
+{
+#if defined(__GNUC__) && !defined(__clang__) && !defined(__CUDA_ARCH__)
+    return std::comp_ellint_1(a);
+#else
+    amrex::ignore_unused(a);
+    AMREX_ALWAYS_ASSERT_WITH_MESSAGE(false, "parser: comp_ellint_1 only supported with gcc and cpu");
+    return 0.0;
+#endif
+}
+
+template <typename T>
+AMREX_GPU_HOST_DEVICE AMREX_NO_INLINE
+T parser_math_comp_ellint_2 (T a)
+{
+#if defined(__GNUC__) && !defined(__clang__) && !defined(__CUDA_ARCH__)
+    return std::comp_ellint_2(a);
+#else
+    amrex::ignore_unused(a);
+    AMREX_ALWAYS_ASSERT_WITH_MESSAGE(false, "parser: comp_ellint_2 only supported with gcc and cpu");
+    return 0.0;
+#endif
+}
+
 template <typename T>
 AMREX_GPU_HOST_DEVICE AMREX_NO_INLINE
 T parser_math_pow (T a, T b) { return std::pow(a,b); }
@@ -297,6 +325,8 @@ parser_call_f1 (enum parser_f1_t type, double a)
     case PARSER_POW_P1:      return a;
     case PARSER_POW_P2:      return a*a;
     case PARSER_POW_P3:      return a*a*a;
+    case PARSER_COMP_ELLINT_1: return parser_math_comp_ellint_1<double>(a);
+    case PARSER_COMP_ELLINT_2: return parser_math_comp_ellint_2<double>(a);
     default:
         amrex::Abort("parser_call_f1: Unknown function ");
         return 0.0;

Src/Base/Parser/AMReX_Parser_Y.cpp

@@ -1209,6 +1209,8 @@ parser_ast_print_f1 (struct parser_f1* f1, std::string const& space, AllPrint& p
     case PARSER_POW_P1:      printer << "POW(,1)\n";     break;
     case PARSER_POW_P2:      printer << "POW(,2)\n";     break;
     case PARSER_POW_P3:      printer << "POW(,3)\n";     break;
+    case PARSER_COMP_ELLINT_1: printer << "COMP_ELLINT_1\n"; break;
+    case PARSER_COMP_ELLINT_2: printer << "COMP_ELLINT_2\n"; break;
     default:
         amrex::AllPrint() << "parser_ast_print_f1: Unknown function " << f1->ftype << "\n";
     }

Src/Base/Parser/amrex_iparser.lex.cpp

@@ -1,9 +1,11 @@
 #ifdef __GNUC__
 #pragma GCC diagnostic ignored "-Wnull-dereference"
 #pragma GCC diagnostic ignored "-Wunreachable-code"
+#pragma GCC diagnostic ignored "-Wsign-compare"
 #elif defined(__clang__)
 #pragma clang diagnostic ignored "-Wnull-dereference"
 #pragma clang diagnostic ignored "-Wunreachable-code"
+#pragma clang diagnostic ignored "-Wsign-compare"
 #endif
 
 #include <amrex_iparser.lex.nolint.H>

Src/Base/Parser/amrex_iparser.lex.h

@@ -281,6 +281,7 @@ typedef int16_t flex_int16_t;
 typedef uint16_t flex_uint16_t;
 typedef int32_t flex_int32_t;
 typedef uint32_t flex_uint32_t;
+typedef uint64_t flex_uint64_t;
 #else
 typedef signed char flex_int8_t;
 typedef short int flex_int16_t;
@@ -360,7 +361,7 @@ typedef struct yy_buffer_state *YY_BUFFER_STATE;
 typedef size_t yy_size_t;
 #endif
 
-extern int yyleng;
+extern yy_size_t yyleng;
 
 extern FILE *yyin, *yyout;
 
@@ -381,7 +382,7 @@ struct yy_buffer_state
 	/* Number of characters read into yy_ch_buf, not including EOB
 	 * characters.
 	 */
-	int yy_n_chars;
+	yy_size_t yy_n_chars;
 
 	/* Whether we "own" the buffer - i.e., we know we created it,
 	 * and can realloc() it to grow it, and should free() it to
@@ -425,7 +426,7 @@ void yypop_buffer_state ( void );
 
 YY_BUFFER_STATE yy_scan_buffer ( char *base, yy_size_t size  );
 YY_BUFFER_STATE yy_scan_string ( const char *yy_str  );
-YY_BUFFER_STATE yy_scan_bytes ( const char *bytes, int len  );
+YY_BUFFER_STATE yy_scan_bytes ( const char *bytes, yy_size_t len  );
 
 void *yyalloc ( yy_size_t  );
 void *yyrealloc ( void *, yy_size_t  );
@@ -482,7 +483,7 @@ FILE *yyget_out ( void );
 
 void yyset_out  ( FILE * _out_str  );
 
-			int yyget_leng ( void );
+			yy_size_t yyget_leng ( void );
 
 char *yyget_text ( void );
 

Src/Base/Parser/amrex_iparser.lex.nolint.H

@@ -303,6 +303,7 @@ typedef int16_t flex_int16_t;
 typedef uint16_t flex_uint16_t;
 typedef int32_t flex_int32_t;
 typedef uint32_t flex_uint32_t;
+typedef uint64_t flex_uint64_t;
 #else
 typedef signed char flex_int8_t;
 typedef short int flex_int16_t;
@@ -411,7 +412,7 @@ typedef struct yy_buffer_state *YY_BUFFER_STATE;
 typedef size_t yy_size_t;
 #endif
 
-extern int yyleng;
+extern yy_size_t yyleng;
 
 extern FILE *yyin, *yyout;
 
@@ -454,7 +455,7 @@ struct yy_buffer_state
 	/* Number of characters read into yy_ch_buf, not including EOB
 	 * characters.
 	 */
-	int yy_n_chars;
+	yy_size_t yy_n_chars;
 
 	/* Whether we "own" the buffer - i.e., we know we created it,
 	 * and can realloc() it to grow it, and should free() it to
@@ -523,8 +524,8 @@ static YY_BUFFER_STATE * yy_buffer_stack = NULL; /**< Stack as an array. */
 
 /* yy_hold_char holds the character lost when yytext is formed. */
 static char yy_hold_char;
-static int yy_n_chars;		/* number of characters read into yy_ch_buf */
-int yyleng;
+static yy_size_t yy_n_chars;		/* number of characters read into yy_ch_buf */
+yy_size_t yyleng;
 
 /* Points to current character in buffer. */
 static char *yy_c_buf_p = NULL;
@@ -551,7 +552,7 @@ static void yy_init_buffer ( YY_BUFFER_STATE b, FILE *file  );
 
 YY_BUFFER_STATE yy_scan_buffer ( char *base, yy_size_t size  );
 YY_BUFFER_STATE yy_scan_string ( const char *yy_str  );
-YY_BUFFER_STATE yy_scan_bytes ( const char *bytes, int len  );
+YY_BUFFER_STATE yy_scan_bytes ( const char *bytes, yy_size_t len  );
 
 void *yyalloc ( yy_size_t  );
 void *yyrealloc ( void *, yy_size_t  );
@@ -607,7 +608,7 @@ static void yynoreturn yy_fatal_error ( const char* msg  );
  */
 #define YY_DO_BEFORE_ACTION \
 	(yytext_ptr) = yy_bp; \
-	yyleng = (int) (yy_cp - yy_bp); \
+	yyleng = (yy_size_t) (yy_cp - yy_bp); \
 	(yy_hold_char) = *yy_cp; \
 	*yy_cp = '\0'; \
 	(yy_c_buf_p) = yy_cp;
@@ -772,7 +773,7 @@ FILE *yyget_out ( void );
 
 void yyset_out  ( FILE * _out_str  );
 
-			int yyget_leng ( void );
+			yy_size_t yyget_leng ( void );
 
 char *yyget_text ( void );
 
@@ -839,7 +840,7 @@ static int input ( void );
 	if ( YY_CURRENT_BUFFER_LVALUE->yy_is_interactive ) \
 		{ \
 		int c = '*'; \
-		int n; \
+		yy_size_t n; \
 		for ( n = 0; n < max_size && \
 			     (c = getc( yyin )) != EOF && c != '\n'; ++n ) \
 			buf[n] = (char) c; \
@@ -1293,7 +1294,7 @@ static int yy_get_next_buffer (void)
 
 	else
 		{
-			int num_to_read =
+			yy_size_t num_to_read =
 			YY_CURRENT_BUFFER_LVALUE->yy_buf_size - number_to_move - 1;
 
 		while ( num_to_read <= 0 )
@@ -1307,7 +1308,7 @@ static int yy_get_next_buffer (void)
 
 			if ( b->yy_is_our_buffer )
 				{
-				int new_size = b->yy_buf_size * 2;
+				yy_size_t new_size = b->yy_buf_size * 2;
 
 				if ( new_size <= 0 )
 					b->yy_buf_size += b->yy_buf_size / 8;
@@ -1365,7 +1366,7 @@ static int yy_get_next_buffer (void)
 
 	if (((yy_n_chars) + number_to_move) > YY_CURRENT_BUFFER_LVALUE->yy_buf_size) {
 		/* Extend the array by 50%, plus the number we really need. */
-		int new_size = (yy_n_chars) + number_to_move + ((yy_n_chars) >> 1);
+		yy_size_t new_size = (yy_n_chars) + number_to_move + ((yy_n_chars) >> 1);
 		YY_CURRENT_BUFFER_LVALUE->yy_ch_buf = (char *) yyrealloc(
 			(void *) YY_CURRENT_BUFFER_LVALUE->yy_ch_buf, (yy_size_t) new_size  );
 		if ( ! YY_CURRENT_BUFFER_LVALUE->yy_ch_buf )
@@ -1468,7 +1469,7 @@ static int yy_get_next_buffer (void)
 
 		else
 			{ /* need more input */
-			int offset = (int) ((yy_c_buf_p) - (yytext_ptr));
+			yy_size_t offset = (yy_c_buf_p) - (yytext_ptr);
 			++(yy_c_buf_p);
 
 			switch ( yy_get_next_buffer(  ) )
@@ -1837,12 +1838,12 @@ YY_BUFFER_STATE yy_scan_string (const char * yystr )
  * 
  * @return the newly allocated buffer state object.
  */
-YY_BUFFER_STATE yy_scan_bytes  (const char * yybytes, int  _yybytes_len )
+YY_BUFFER_STATE yy_scan_bytes  (const char * yybytes, yy_size_t  _yybytes_len )
 {
 	YY_BUFFER_STATE b;
 	char *buf;
 	yy_size_t n;
-	int i;
+	yy_size_t i;
 
 	/* Get memory for full buffer, including space for trailing EOB's. */
 	n = (yy_size_t) (_yybytes_len + 2);
@@ -1884,7 +1885,7 @@ static void yynoreturn yy_fatal_error (const char* msg )
 	do \
 		{ \
 		/* Undo effects of setting up yytext. */ \
-        int yyless_macro_arg = (n); \
+        yy_size_t yyless_macro_arg = (n); \
         YY_LESS_LINENO(yyless_macro_arg);\
 		yytext[yyleng] = (yy_hold_char); \
 		(yy_c_buf_p) = yytext + yyless_macro_arg; \
@@ -1924,7 +1925,7 @@ FILE *yyget_out  (void)
 /** Get the length of the current token.
  * 
  */
-int yyget_leng  (void)
+yy_size_t yyget_leng  (void)
 {
         return yyleng;
 }

Src/Base/Parser/amrex_iparser.tab.h

@@ -1,4 +1,4 @@
-/* A Bison parser, made by GNU Bison 3.7.5.  */
+/* A Bison parser, made by GNU Bison 3.8.2.  */
 
 /* Bison interface for Yacc-like parsers in C
 
@@ -16,7 +16,7 @@
    GNU General Public License for more details.
 
    You should have received a copy of the GNU General Public License
-   along with this program.  If not, see <http://www.gnu.org/licenses/>.  */
+   along with this program.  If not, see <https://www.gnu.org/licenses/>.  */
 
 /* As a special exception, you may create a larger work that contains
    part or all of the Bison parser skeleton and distribute that work
@@ -105,6 +105,8 @@ typedef union AMREX_IPARSERSTYPE AMREX_IPARSERSTYPE;
 
 extern AMREX_IPARSERSTYPE amrex_iparserlval;
 
+
 int amrex_iparserparse (void);
 
+
 #endif /* !YY_AMREX_IPARSER_AMREX_IPARSER_TAB_H_INCLUDED  */