R devel - Aug 2011 - How to safely using OpenMP pragma inside a .C() function?

I am trying to parallelize part of a C function that is called from R (via
.C) using OpenMP's "parallel for" pragma. I get mixed results:
some runs
finish with no problem, but some lead to R crashing after issuing a long
error message involving memory violations.

I found this post, which describes how a .Call() function can be made to
avoid crashing R by raising the stack limit:

http://www.r-bloggers.com/using-openmp-ized-c-code-with-r/

However, trying this in my .C function doesn't help things. Any
suggestions/tips on whether I can safely use OpenMP inside a .C function,
and if yes how?

Thank you,
Alireza Mahani

--
View this message in context:
http://r.789695.n4.nabble.com/How-to-safely-using-OpenMP-pragma-inside-a-C-function-tp3777036p3777036.html
Sent from the R devel mailing list archive at Nabble.com.

On Aug 29, 2011, at 7:48 PM, Alireza Mahani wrote:
> I am trying to parallelize part of a C function that is called from R (via
> .C) using OpenMP's "parallel for" pragma. I get mixed
results: some runs
> finish with no problem, but some lead to R crashing after issuing a long
> error message involving memory violations.
> 
You'll have to provide the code. In general it works (even R uses it
itself), but there are strict requirements (no R API calls) that you must adhere
to.

> I found this post, which describes how a .Call() function can be made to
> avoid crashing R by raising the stack limit:
> 
> http://www.r-bloggers.com/using-openmp-ized-c-code-with-r/
> 
I skimmed through the post and all of the examples are broken - they will only
work (incidentally) as R internals, not officially (and they are unnecessary
inefficient).

> However, trying this in my .C function doesn't help things. Any
> suggestions/tips on whether I can safely use OpenMP inside a .C function,
> and if yes how?
> 
There are issues with OpenMP on some platforms in general (in fact pretty much
every platform had some issue at some point in time), but apart from those you
only have to make sure that you declare shared/private variables properly and
don't use *any* R API calls in the parallel part (this includes things like
LENGTH, REAL, ...).

Cheers,
Simon


> Thank you,
> Alireza Mahani
> 
> --
> View this message in context:
http://r.789695.n4.nabble.com/How-to-safely-using-OpenMP-pragma-inside-a-C-function-tp3777036p3777036.html
> Sent from the R devel mailing list archive at Nabble.com.
> 
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
> 
>

I have no R API calls inside the parallelized block. I will work on creating
a self-contained example and post it for your review. Thanks! -Alireza

--
View this message in context:
http://r.789695.n4.nabble.com/How-to-safely-use-OpenMP-pragma-inside-a-C-function-tp3777036p3778482.html
Sent from the R devel mailing list archive at Nabble.com.

Simon,

I found that files R-2.13.1/src/library/stats/src/distance.c and
R-2.13.1/src/main/array.c have openmp code (example below). I have couple
questions regarding best practices when using R internals and openmp. 

Can we use R-2.13.1/src/library/stats/src/distance.c and
R-2.13.1/src/main/array.c as an example how to interact with R code and R
internals?
What are my other options if I want to work with SEXP structures in my
parallel code?

Thank you
Regards

============
#ifdef HAVE_OPENMP
        /* This gives a spurious -Wunused-but-set-variable error */
        if (R_num_math_threads > 0) 
            nthreads = R_num_math_threads;
        else 
            nthreads = 1; /* for now */
#pragma omp parallel for num_threads(nthreads) default(none) \
    private(j, i, ix, rx) \
    firstprivate(x, ans, n, p, type, cnt, sum, \
                 NaRm, keepNA, R_NaReal, R_NaInt, OP)
#endif
        for (j = 0; j < p; j++) {
            switch (type) {
            case REALSXP:
                rx = REAL(x) + n*j; 
                if (keepNA)
                    for (sum = 0., i = 0; i < n; i++) sum += *rx++;
                else {
                    for (cnt = 0, sum = 0., i = 0; i < n; i++, rx++)
                        if (!ISNAN(*rx)) {cnt++; sum += *rx;}
                        else if (keepNA) {sum = NA_REAL; break;}
                }    
                break;
            case INTSXP:
                ix = INTEGER(x) + n*j; 
                for (cnt = 0, sum = 0., i = 0; i < n; i++, ix++)
                    if (*ix != NA_INTEGER) {cnt++; sum += *ix;}
                    else if (keepNA) {sum = NA_REAL; break;}
                break;
            case LGLSXP:
                ix = LOGICAL(x) + n*j; 
                for (cnt = 0, sum = 0., i = 0; i < n; i++, ix++)
                    if (*ix != NA_LOGICAL) {cnt++; sum += *ix;}
                    else if (keepNA) {sum = NA_REAL; break;}
                break;
            default:
                /* we checked the type above, but be sure */
                UNIMPLEMENTED_TYPEt("do_colsum", type);
            }
            if (OP == 1) {
                if (cnt > 0) sum /= cnt; else sum = NA_REAL;
            }
            REAL(ans)[j] = sum;
        }


--
View this message in context:
http://r.789695.n4.nabble.com/How-to-safely-use-OpenMP-pragma-inside-a-C-function-tp3777036p3779214.html
Sent from the R devel mailing list archive at Nabble.com.