##### https://github.com/cran/nacopula
Tip revision: e24ab1e
rFrank.c
``````/*
Copyright (C) 2010 Marius Hofert and Martin Maechler

This program is free software; you can redistribute it and/or modify it under
Foundation; either version 3 of the License, or (at your option) any later
version.

This program is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
FOR A PARTICULAR PURPOSE. See the 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/>.
*/

#include <Rmath.h>

#include "nacopula.h"

/**
* Sample V ~ F with Laplace-Stieltjes transform
* (1-(1-exp(-t)*(1-e^(-theta1)))^alpha)/(1-e^(-theta0))
* via the algorithm of Hofert (2010). C version.
* @param p parameter 1-e^(-theta1)
* @param alpha parameter theta0/theta1 in (0,1]
* @param iAlpha 1-alpha
* @param theta0_le_1 in {0,1} with 1 if and only if theta0 <= 1
* @return a random variate from F
* @author Marius Hofert, Martin Maechler
*/
double rFFrank(double p, double alpha, double iAlpha /**< == 1 - alpha */,
int theta0_le_1){
double U, V;
if(theta0_le_1) {
do {
U = unif_rand();
V = rLog(p);
} while (U*(V-alpha) > 1./beta(V, iAlpha));

} else {
double gamma_1_a = gammafn(iAlpha);
do {
U = unif_rand();
V = rFJoe(alpha, iAlpha, gamma_1_a);
} while(U > pow(p, V-1.));
}
return V;
}

/**
* Vectorize rFFrank. Generate a vector of variates
* V ~ F with Laplace-Stieltjes transform
* (1-(1-exp(-t)*(1-e^(-theta1)))^alpha)/(1-e^(-theta0)).
* @param V vector of random variates from F (result)
* @param n length of the vector V
* @param theta_0 parameter theta0 in (0,infinity)
* @param theta_1 parameter theta1 in [theta0, infinity)
* @return none
* @author Marius Hofert, Martin Maechler
*/
void rFFrank_vec(double *V, const int n, const double theta_0,
const double theta_1){
double p  =  - expm1(-theta_1),
alpha = theta_0 / theta_1, iAlpha = (theta_1 - theta_0) / theta_1;
int th_0_le_1 = (theta_0 <= theta_1);

if(n >= 1) {
GetRNGstate();

for(int i=0; i < n; i++)
V[i] = rFFrank(p, alpha, iAlpha, th_0_le_1);

PutRNGstate();
}
return;
}

/**
* Generate a vector of variates V ~ F with Laplace-Stieltjes transform
* (1-(1-exp(-t)*(1-e^(-theta1)))^alpha)/(1-e^(-theta0)).
* @param n sample size
* @param theta_0 parameter theta0 in (0,infinity)
* @param theta_1 parameter theta1 in [theta0, infinity)
* @return vector of random variates V
* @author Martin Maechler
*/
SEXP rFFrank_c(SEXP n_, SEXP theta_0_, SEXP theta_1_){
int n = asInteger(n_);
double theta_0 = asReal(theta_0_),
theta_1 = asReal(theta_1_);
SEXP res = PROTECT(allocVector(REALSXP, n));
if(n >= 1)
rFFrank_vec(REAL(res), n, theta_0, theta_1);
UNPROTECT(1);
return res;
}
``````