A simple CLI arbitrary-precision Mandelbrot set calculator (for a single point)

Question

I wrote the following C code that takes a couple of arguments on the command line and calculates a given number of iterations of the Mandelbrot set complex mapping using MPFR arbitrary precision arithmetic. Is it well-written? Does it have undefined behaviour or is unsafe?

Here is the code:

Makefile

CC := gcc
O := 1
CFLAGS := -Wall -Wpedantic -std=gnu17 -O$(O) $(XCFLAGS)
XCFLAGS :=

mpfrtest: mpfrtest.o
    $(CC) $(CFLAGS) -lm -lmpfr $^ -o $@

mpfrtest.o: mpfrtest.c
    $(CC) $(CFLAGS) -c $^ -o $@

clean:
    rm -f *.o
    rm -f mpfrtest

mpfrtest.c

#include <complex.h>
#include <errno.h>
#include <inttypes.h>
#include <math.h>
#include <mpfr.h>
#include <stdio.h>
#include <stdlib.h>

#define ESCAPE_NORM 4
#define ROUND MPFR_RNDN

size_t orbit(const mpfr_t c_re, const mpfr_t c_im, size_t M,
             double complex *output, mpfr_prec_t prec, mpfr_rnd_t round,
             void (*callback)(size_t)) {
    mpfr_t z_re, z_im;
    mpfr_inits2(prec, z_re, z_im, NULL);
    mpfr_set_ui(z_re, 0, round);
    mpfr_set_ui(z_im, 0, round);
    mpfr_t re_tmp, im_tmp;
    mpfr_inits2(prec, re_tmp, im_tmp, NULL);
    double interval_d = ceil(sqrt((double)M));
    size_t interval;
    if (interval_d < 1) {
        interval = 1;
    } else {
        interval = (size_t)interval_d;
    }
    output[0] = 0;
    size_t m = 0;
    for (; m < M && mpfr_cmp_ui(re_tmp, ESCAPE_NORM) <= 0; m++) {
        mpfr_sqr(re_tmp, z_re, round);
        mpfr_sqr(im_tmp, z_im, round);
        mpfr_sub(im_tmp, im_tmp, c_re, round);
        mpfr_sub(re_tmp, re_tmp, im_tmp, round);
        mpfr_mul(z_im, z_re, z_im, round);
        mpfr_mul_ui(z_im, z_im, 2, round);
        mpfr_add(z_im, z_im, c_im, round);
        mpfr_set(z_re, re_tmp, round);
        output[m] = mpfr_get_d(z_re, round) + mpfr_get_d(z_im, round) * I;
        mpfr_sqr(re_tmp, z_re, round);
        mpfr_sqr(im_tmp, z_im, round);
        mpfr_add(re_tmp, re_tmp, im_tmp, round);
        if (m % interval == 0) {
            (*callback)(m);
        }
    }
    mpfr_clears(re_tmp, im_tmp, z_re, z_im, NULL);
    return m;
}

int digitsiz(size_t value) {
    value++;
    int count = 0;
    while (value >= 10) {
        value /= 10;
        count++;
    }
    return count;
}

void progress(size_t m) { fprintf(stderr, "computed iteration %zu\n", m); }

int main(int argc, char **argv) {
    if (argc != 5) {
        fprintf(stderr, "wrong argument count %d, should be 4\n", argc - 1);
        fprintf(
            stderr,
            "usage: %s <Re[c]> <Im[c]> <precision> <maximum iteration count>\n",
            argv[0]);
        exit(-5);
    }
    char *endptr;
    errno = 0;
    signed long prec_l = strtol(argv[3], &endptr, 10);
    if (argv[3][0] == '\0' || *endptr != '\0') {
        fprintf(stderr, "failed to parse argument precision as a number\n");
        exit(-3);
    }
    if (errno == ERANGE) {
        fprintf(stderr, "argument precision is too large\n");
        exit(-3);
    }
    if (prec_l <= 0) {
        fprintf(stderr, "argument precision must be positive\n");
        exit(-3);
    }
    mpfr_prec_t prec = (mpfr_prec_t)prec_l;
    errno = 0;
    uintmax_t iter_umax = strtoumax(argv[4], &endptr, 10);
    if (argv[4][0] == '\0' || *endptr != '\0') {
        fprintf(
            stderr,
            "failed to parse argument maximum iteration count as a number\n");
        exit(-4);
    }
    if (iter_umax > (uintmax_t)SIZE_MAX || errno == ERANGE) {
        fprintf(stderr, "argument maximum iteration count is too large\n");
        exit(-4);
    }
    size_t iter = (size_t)iter_umax;
    mpfr_t c_re;
    mpfr_t c_im;
    mpfr_inits2(prec, c_re, c_im, NULL);
    if (mpfr_set_str(c_re, argv[1], 10, ROUND) != 0) {
        fprintf(stderr, "failed to parse argument Re[c] as a number\n");
        mpfr_clears(c_re, c_im, NULL);
        exit(-4);
    };
    if (mpfr_set_str(c_im, argv[2], 10, ROUND) != 0) {
        fprintf(stderr, "failed to parse argument Im[c] as a number\n");
        mpfr_clears(c_re, c_im, NULL);
        exit(-5);
    }
    double complex *result = malloc(iter * sizeof(double complex));
    if (result == NULL) {
        perror("orbit array allocation");
        mpfr_clears(c_re, c_im, NULL);
        exit(-6);
    }
    size_t m = orbit(c_re, c_im, iter, result, prec, ROUND, progress);
    if (m == SIZE_MAX) {
        printf("iterations = max");
    } else {
        printf("iterations = %zu\n", m);
    }
    int width = digitsiz(m);
    for (size_t i = 0; i < m; i++) {
        printf("orbit[%0*zu] = %13.6e + %13.6ei\n", width, i, creal(result[i]),
               cimag(result[i]));
    }
    mpfr_clears(c_re, c_im, NULL);
    free(result);
}

Answer 1

Only for the pedantic

strtou...() functions have a quirk in that strtoumax("-1", &endptr, 10) is some large positive value without setting errno.

If that is OK, 'nuff said.

If code wants to detect negative strings, consider strtomax() and test the range.

// Demonstrative unchecked code.

// Like strtoumax()`, yet caps negative conversions rather than wrap.
// This does allow "-0" to quietly return 0 with no setting of `errno`.
uintmax_t strtoumax_alt(const char * restrict nptr, char ** restrict endptr, int base) {
  // Usually `errno = 0` set by calling code.
  intmax_t i = strtomax(nptr, endptr, base);
  if (i > INTMAX_MAX - 1) {
    if (errno == ERANGE) {
      errno = 0;
    }
    return strtoumax(nptr, endptr, base);
  }
  if (i < 0) {
    errno = ERANGE;
    return 0;
  }
  return (uintmax_t) i;
}

errno complicates things.^*1

---

argc == 0 possible

Following code is then UB.

if (argc != 5) {
    ...
    fprintf(stderr, "usage: %s ...\n", argv[0]);

---

^*1

strto...() functions are most clearly used if the calling code first performs errno = 0;

This strtoumax_alt(), with its two step of approach of strtointmax()/strtouintmax() like-wise works best if the caller first performs errno = 0;

Yet if we want to better handle cases where strtoumax_alt() is called with a non-zero errno and regular strtoumax() does not set errno with the given parameters, then consider restoring errno before the 2nd strto...() function. So if the calling code is somehow accumulating errors in errno, this restoration will better support that.

uintmax_t strtoumax_alt(const char * restrict nptr, char ** restrict endptr, int base) {
  int errno_original = errno;
  intmax_t i = strtomax(nptr, endptr, base);
  if (i > INTMAX_MAX - 1) {
    errno = errno_original;
    return strtoumax(nptr, endptr, base);
  }
  ...

I suppose now there may be a race condition if 2 threads execute strtoumax_alt(), yet often such threads roll their own copy of errno.

Messing with global objects for error indication is tricky.

How I wish C afforded a direct return value and error flag feature (without a struct).

Answer 2

General Observations

The good, the call to malloc() is tested for success, the code is properly indented and there are symbolic constants rather than raw numbers in most cases. Symbols that coude be added are EXPECTED_ARG_COUNT 5.

The code is too tightly packed, it needs more vertical spacing. It is very difficult to read and follow. There should be blank lines in orbit() and main() to separate the logic blocks (actually in all the functions but it is most notable in orbit() and main()).

The code is not user friendly, it expects all input to be in a specific order. It would be easier for users if each expected argument was preceded by a flag (for example --precision <precision value>).

Since the code is generating Mandelbrots I'm not sure why it isn't tied to graphic output.

Use the Compiler to Find Possible Issues in the Code

I always compile with these flags to find possible issues in the code -Wall -Wextra -pedantic.

In this case the gcc compiler generates the following warning messages:

main.c: In function 'progress':
main.c:62:64: warning: unknown conversion type character 'z' in format [-Wformat=]
void progress(size_t m) { fprintf(stderr, "computed iteration %zu\n", m); }
main.c:62:43: warning: too many arguments for format [-Wformat-extra-args]
void progress(size_t m) { fprintf(stderr, "computed iteration %zu\n", m); }
main.c:62:64: warning: unknown conversion type character 'z' in format [-Wformat=]
void progress(size_t m) { fprintf(stderr, "computed iteration %zu\n", m); }
main.c:62:43: warning: too many arguments for format [-Wformat-extra-args]
void progress(size_t m) { fprintf(stderr, "computed iteration %zu\n", m); }
main.c: In function 'main':
main.c:126:31: warning: unknown conversion type character 'z' in format [-Wformat=]
printf("iterations = %zu\n", m);
main.c:126:16: warning: too many arguments for format [-Wformat-extra-args]
printf("iterations = %zu\n", m);
main.c:126:31: warning: unknown conversion type character 'z' in format [-Wformat=]
printf("iterations = %zu\n", m);
main.c:126:16: warning: too many arguments for format [-Wformat-extra-args]
printf("iterations = %zu\n", m);
main.c:130:26: warning: unknown conversion type character 'z' in format [-Wformat=]
printf("orbit[%0zu] = %13.6e + %13.6ei\n", width, i, creal(result[i]),
main.c:130:37: warning: format '%e' expects argument of type 'double', but argument 3 has type 'size_t' {aka 'long long unsigned int'} [-Wformat=]
printf("orbit[%0zu] = %13.6e + %13.6ei\n", width, i, creal(result[i]),
%13.6I64d
main.c:130:16: warning: too many arguments for format [-Wformat-extra-args]
printf("orbit[%0zu] = %13.6e + %13.6ei\n", width, i, creal(result[i]),
main.c:130:26: warning: unknown conversion type character 'z' in format [-Wformat=]
printf("orbit[%0zu] = %13.6e + %13.6ei\n", width, i, creal(result[i]),
main.c:130:37: warning: format '%e' expects argument of type 'double', but argument 3 has type 'size_t' {aka 'long long unsigned int'} [-Wformat=]
printf("orbit[%0zu] = %13.6e + %13.6ei\n", width, i, creal(result[i]),
%13.6I64d
main.c:130:16: warning: too many arguments for format [-Wformat-extra-args]
printf("orbit[%0zu] = %13.6e + %13.6ei\n", width, i, creal(result[i]),

Don't Call exit(int status) From main()

Since the main() function is the entry point for the program there is no need to call the exit() function, all that is needed is to return with the program status.

Example
Instead of

        exit(-3);

the code should be

       return(-3);

System Provided Exit Values

All of the exit() calls contain raw numbers, there are 2 system provided exit values in the C programming language (also the C++ programming language). These are EXIT_SUCCESS and EXIT_FAILURE. These values should be used instead of raw numbers since they are system defined and system dependent.

Complexity

The function main() is too complex (does too much). As programs grow in size the use of main() should be limited to calling functions that parse the command line, calling functions that set up for processing, calling functions that execute the desired function of the program, and calling functions to clean up after the main portion of the program.

There is also a programming principle called the Single Responsibility Principle that applies here. The Single Responsibility Principle states:

that every module, class, or function should have responsibility over a single part of the functionality provided by the software, and that responsibility should be entirely encapsulated by that module, class or function.

There are at least 3 possible functions in main().

• Parse the command line
• Set up the variables for processing
• Call orbit()
• Output the results