Hello,
My name is Ick Hoon Jin and I am Ph. D. student in Texas A & M Univ..
When I run the C embedded in R in the Linux system, I confront the
following error after 6,000 iteration. By googling I found this error is
from the problem in C.

*** caught segfault ***
address (nil), cause 'memory not mapped'

My C code is following:

---------------------------------------------------------------------------
#include<stdio.h>
#include<stdlib.h>
#include<R.h>
#include<Rmath.h>
#include<malloc.h>

#define Min(X,Y) ((X) < (Y) ? (X) : (Y))

void SecondMH(int *n, int *X, int *length_X, double *theta_new, int
*n_para, double *out);
int MH_Result(double *theta_new, int n_para, int *M1, int *M2, int n);
double Count_Edges(int *X, int n);
double Cal_GWD(int *X, int n);
int Cal_Degree(int *M3, int n);
double Inner_Product( double *vector1, double *vector2, int length);
int rdtsc();

void SecondMH(int *n, int *X, int *length_X, double *theta_new, int
*n_para, double *out)
{
        int i, j, a, k;
        int *M1, *M2;
        M1 = (int *)calloc(*length_X, sizeof(int));
        M2 = (int *)calloc(*length_X, sizeof(int));

        for(i = 1; i < *n; i++)
        {
                for(j = 0; j < i; j++)
                {
                        for(k = 0; k < *length_X; k++ )
                        {
                                M1[k] = X[k];
                                M2[k] = X[k];
                        }

                        if(X[i * *n + j] == 1)
                        {
                                M1[i * *n + j] = 0;
                                a = MH_Result(theta_new, *n_para, M1, M2, *n);
                                if(a == 1)
                                {
                                        X[i * *n + j] = 0;
                                }
                        }
                        else
                        {
                                M1[i * *n + j] = 1;
                                a = MH_Result(theta_new, *n_para, M1, M2, *n);
                                if(a == 1)
                                {
                                        X[i * *n + j] = 1;
                                }
                        }
                }
        }

        for(i = 1; i < *n; i++)
        {
                for(j = 0; j < i; j++)
                {
                        X[j * *n + i] = 0;
                }
        }

        for(i = 0; i < *length_X; i++)
        {
                out[i] = (double)X[i];
        }

        free(M1);
        free(M2);

        return;
}

int MH_Result(double *theta_new, int n_para, int *M1, int *M2, int n)
{
        double *M1_STAT, *M2_STAT;
        double pi_Num, pi_Denom, MH_Ratio, v;

        M1_STAT = (double *)calloc( n_para, sizeof( double ) );
        M2_STAT = (double *)calloc( n_para, sizeof( double ) );

        M1_STAT[0] = Count_Edges(M1, n);
        M2_STAT[0] = Count_Edges(M2, n);

        M1_STAT[1] = Cal_GWD(M1, n);
        M2_STAT[1] = Cal_GWD(M2, n);

        pi_Num = Inner_Product(theta_new, M1_STAT, n_para);
        pi_Denom = Inner_Product(theta_new, M2_STAT, n_para);
        MH_Ratio = pi_Num - pi_Denom;

        srand(rdtsc());
        v = (double)rand() / ( (double)RAND_MAX + (double)1 );

        if( log( v ) < Min( 0, MH_Ratio ) )
                return 1;
        else
                return 0;

        free(M1_STAT);
        free(M2_STAT);
}

double Count_Edges(int *X, int n)
{
        double temp;
        int i, j;

        temp = 0;

        for(i = 1; i < n; i++)
        {
                for(j = 0; j < i; j++)
                {
                        temp += X[i * n + j];
                }
        }

        return temp;
}

double Cal_GWD(int *X, int n)
{

        int *M3, *Degree, *D_Y;
        int i, j;
        double theta = 0.25;
        double GWD = 0.00;

        M3 = (int *)calloc(n, sizeof(int));
        Degree = (int *)calloc(n, sizeof(int));
        D_Y = (int *)calloc((n - 1), sizeof(int));

        for(i = 0; i < n; i++)
        {
                for(j = 0; j < n; j++)
                {
                        M3[j] = X[i * n + j];
                }

                Degree[i] = Cal_Degree(M3, n);
        }

        for(i = 0; i < n - 1; i++)
        {
                D_Y[i] = 0;

                for(j = 0; j < n; j++)
                {
                        if(Degree[j] = i + 1)
                        {
                                D_Y[i] += 1;
                        }
                }
        }

        for(i = 0; i < n-1; i++)
        {
                        GWD += ( 1 - pow( (1 - exp(theta*(-1))), (i + 1) ) ) * 
(double)D_Y[i];
        }

        GWD = GWD * exp(theta);

        free(M3);
        free(Degree);
        free(D_Y);

        return GWD;
}

int Cal_Degree(int *M3, int n)
{
        int i, result;

        result = 0;

        for(i = 0; i < n; i++)
        {
                result += M3[i];
        }

        return result;
}

double Inner_Product(double *vector1, double *vector2, int length)
{
        int i;
        double result;

        result = 0;
        for(i = 0; i < length; i++)
        {
                result += vector1[i] * vector2[i];
        }

        return result;
}

int rdtsc()
{
    __asm__ __volatile__("rdtsc");
}
-------------------------------------------------------------------------

and my R-code to call C function is following;

#-------------------------------------------------------------------------------

# Function for Second Gibbs with Metropolis-Hasting Step

#-------------------------------------------------------------------------------

dyn.load("SecondMH.so")

second.MH<-function(Adj.Matrix,theta.new,nodes){

  n<-as.integer(nodes)

  Y<-as.matrix.network(X,matrix.type="adjacency",directed=FALSE)

  Z<-as.integer(as.vector(Y))

  length.Z<-as.integer(length(Z))

  theta.star<-as.double(as.vector(theta.new))

  n.para<-as.integer(length(theta.star))

  
result<-.C("SecondMH",n,Z,length.Z,theta.star,n.para,out=as.double(rep(0.0,length.Z)))

  X.new<-matrix(result$out,n,n)
  Adj.Matrix<-network(X.new, directed=FALSE)

  return(Adj.Matrix)}


--------------------------------------------------------------------------

Until this time, I have used only R and I am really a beginner for C
programming. How can I solve the segfault error? Please answer it.
Thank you very much.
Sincerely,



Jin, Ick Hoon
Ph. D. Student, Department of Statistics
Texas A&M University, College Station, TX

______________________________________________
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel

Reply via email to