#Use what C++ provides If you use the C++ versions of the headers:

#include <cstdlib>
#include <cstdio> // or better, <iostream>
#include <cstring>

You'll get the names safely in the std namespace. I see you include <cmath> but assume its names are exported to the global namespace - that's an error. You can't depend on that happening (it's just an artefact of your implementation).

Also, prefer new/delete and new[]/delete[] over malloc()/free(). You might be able to avoid memory management completely by using the standard containers; that's a Good Thing.

std::max() and std::min() are safer choices compared to your maxd() and mind() macros.

#Avoid I/O in parallel code These prints may get interleaved in unpredictable ways:

#pragma omp parallel for
for (int i = 0; i < cur_root; i += 2)
  {
    fprintf (temp, "%lf %lf %lf\n",
             roots[i], roots[i+1],
             calculateDensityScalarFromDensity(min, max, densities[i / 2])); 
  }

You might do something like

#pragma omp parallel for
for (int i = 0; i < cur_root; i += 2)
  {
    auto density = calculateDensityScalarFromDensity(min, max, densities[i / 2])
#pragma omp critical
    fprintf (temp, "%lf %lf %lf\n",
             roots[i], roots[i+1],
             density); 
  }

#Use what C++ provides If you use the C++ versions of the headers:

#include <cstdlib>
#include <cstdio> // or better, <iostream>
#include <cstring>

You'll get the names safely in the std namespace. I see you include <cmath> but assume its names are exported to the global namespace - that's an error. You can't depend on that happening (it's just an artefact of your implementation).

Also, prefer new/delete and new[]/delete[] over malloc()/free(). You might be able to avoid memory management completely by using the standard containers; that's a Good Thing.

std::max() and std::min() are safer choices compared to your maxd() and mind() macros.

#Avoid I/O in parallel code These prints may get interleaved in unpredictable ways:

#pragma omp parallel for
for (int i = 0; i < cur_root; i += 2)
  {
    fprintf (temp, "%lf %lf %lf\n",
             roots[i], roots[i+1],
             calculateDensityScalarFromDensity(min, max, densities[i / 2])); 
  }

You might do something like

#pragma omp parallel for
for (int i = 0; i < cur_root; i += 2)
  {
    auto density = calculateDensityScalarFromDensity(min, max, densities[i / 2])
#pragma omp critical
    fprintf (temp, "%lf %lf %lf\n",
             roots[i], roots[i+1],
             density); 
  }

#Bin your results I don't mean to throw everything away! The output is a count of results at each pixel location (possibly with some small blurring? I'm not quite sure). Having obtained values, the only use for them is to increment the count of the corresponding bin - it's like obtaining a histogram. That reduces your storage requirements, and your output is almost formatted for you.

#Use symmetry As noted in comments, we can get away with recording only the positive roots, and produce the rest of the output by reflection.