SEXP sampler_glue_C_dist(
    SEXP sampler_name, SEXP sampler_context, SEXP log_dens_name,
    SEXP dist_context, SEXP x0, SEXP sample_size, SEXP tuning) {

  // Locate symbol for sampler function.

  const char *sampler_str = CHAR(STRING_ELT(sampler_name,0));
  sampler_t *sampler_fp = (sampler_t*)R_FindSymbol(sampler_str, "", NULL);
  if (sampler_fp==NULL)
    error("Cannot locate symbol \"%s\".", sampler_str);

  // Locate symbol for log density.

  const char *log_dens_str = CHAR(STRING_ELT(log_dens_name,0));
  log_density_t *log_dens_fp =
    (log_density_t*)R_FindSymbol(log_dens_str, "", NULL);
  if (log_dens_fp==NULL)
    error("Cannot locate symbol \"%s\".", log_dens_str);

  // Define a stub function to keep track of the number of function calls.

  int ndim = length(x0);
  C_stub_context_t stub_context =
    { .ds = { .log_dens=log_dens_fp, .ndim=ndim, .context=dist_context },
      .evals=0, .grads=0 };
  SEXP raw_context;
  PROTECT(raw_context = void_as_raw(&stub_context));
  dist_t stub_ds = { .log_dens=C_log_density_stub_func,
                     .context=raw_context, .ndim=ndim };

  // Create a matrix to store the states in and call the sampler.

  SEXP X;
  PROTECT(X = allocMatrix(REALSXP, *REAL(sample_size), ndim));
  GetRNGstate();
  sampler_fp(sampler_context, &stub_ds, REAL(x0), *REAL(sample_size),
             *REAL(tuning), REAL(X));
  PutRNGstate();

  // Construct the result to return.

  const char *result_names[] = { "X", "evals", "grads", "" };
  SEXP result;
  PROTECT(result = mkNamed(VECSXP, result_names));
  SET_VECTOR_ELT(result, 0, X);
  SET_VECTOR_ELT(result, 1, ScalarInteger(stub_context.evals));
  SET_VECTOR_ELT(result, 2, ScalarInteger(stub_context.grads));
  UNPROTECT(3);
  return result;
}

/*
  This is the routine associated with the getNativeSymbolInfo()
  function and it takes the name of a symbol and optionally an
  object identifier (package usually) in which to restrict the search
  for this symbol. It resolves the symbol and returns it to the caller
  giving the symbol address, the package information (i.e. name and
  fully qualified shared object name). If the symbol was explicitly
  registered (rather than dynamically resolved by R), then we pass
  back that information also, giving the number of arguments it
  expects and the interface by which it should be called.
  The returned object has class NativeSymbol. If the symbol was
  registered, we add a class identifying the interface type
  for which it is intended (i.e. .C(), .Call(), etc.)
 */
SEXP attribute_hidden
R_getSymbolInfo(SEXP sname, SEXP spackage, SEXP withRegistrationInfo)
{
    const void *vmax = vmaxget();
    const char *package, *name;
    R_RegisteredNativeSymbol symbol = {R_ANY_SYM, {NULL}, NULL};
    SEXP sym = R_NilValue;
    DL_FUNC f = NULL;

    package = "";

    name = translateChar(STRING_ELT(sname, 0));

    if(length(spackage)) {
	if(TYPEOF(spackage) == STRSXP)
	    package = translateChar(STRING_ELT(spackage, 0));
	else if(TYPEOF(spackage) == EXTPTRSXP &&
		R_ExternalPtrTag(spackage) == install("DLLInfo")) {
	    f = R_dlsym((DllInfo *) R_ExternalPtrAddr(spackage), name, &symbol);
	    package = NULL;
	} else
	    error(_("must pass package name or DllInfo reference"));
    }

    if(package)
	f = R_FindSymbol(name, package, &symbol);

    if(f)
	sym = createRSymbolObject(sname, f, &symbol,
				  LOGICAL(withRegistrationInfo)[0]);

    vmaxset(vmax);
    return sym;
}

SEXP attribute_hidden
do_getSymbolInfo(SEXP call, SEXP op, SEXP args, SEXP env)
{
    const char *package = "", *name;
    R_RegisteredNativeSymbol symbol = {R_ANY_SYM, {NULL}, NULL};
    SEXP sym = R_NilValue;
    DL_FUNC f = NULL;

    checkArity(op, args);
    SEXP sname = CAR(args), spackage = CADR(args), 
	withRegistrationInfo = CADDR(args);

    name = translateChar(STRING_ELT(sname, 0));
    if(length(spackage)) {
	if(TYPEOF(spackage) == STRSXP)
	    package = translateChar(STRING_ELT(spackage, 0));
	else if(TYPEOF(spackage) == EXTPTRSXP &&
		R_ExternalPtrTag(spackage) == install("DLLInfo")) {
	    f = R_dlsym((DllInfo *) R_ExternalPtrAddr(spackage), name, &symbol);
	    package = NULL;
	} else
	    error(_("must pass package name or DllInfo reference"));
    }
    if(package)
	f = R_FindSymbol(name, package, &symbol);
    if(f)
	sym = createRSymbolObject(sname, f, &symbol,
				  LOGICAL(withRegistrationInfo)[0]);
    return sym;
}

QuartzFunctions_t *getQuartzFunctions(void) {
    if (qfn) return qfn;
    {
	QuartzFunctions_t *(*fn)(void);
	fn = (QuartzFunctions_t *(*)(void)) R_FindSymbol("getQuartzAPI", "grDevices", NULL);
	if (!fn) {
	    /* we need to load grDevices - not sure if this is the best way, though ... */
	    SEXP call = lang2(install("library"), install("grDevices"));
	    PROTECT(call);
	    eval(call, R_GlobalEnv);
	    UNPROTECT(1);
	    fn = (QuartzFunctions_t *(*)(void)) R_FindSymbol("getQuartzAPI", "grDevices", NULL);
	    if (!fn) error(_("unable to load Quartz"));
	}
	return fn();
    }
}

static void R_de_Init(void)
{
    static int de_init = 0;

    if(de_init > 0) return;
    if(de_init < 0) error(_("X11 dataentry cannot be loaded"));

    de_init = -1;
    if(strcmp(R_GUIType, "none") == 0) {
	warning(_("X11 is not available"));
	return;
    }
    int res = R_moduleCdynload("R_de", 1, 1);
    if(!res) error(_("X11 dataentry cannot be loaded"));
    de_ptr->de = (R_X11DataEntryRoutine) 
	R_FindSymbol("in_RX11_dataentry", "R_de", NULL);
    de_ptr->dv = (R_X11DataViewer) 
	R_FindSymbol("in_R_X11_dataviewer", "R_de", NULL);
    de_init = 1;
    return;
}

/*********************
 void WtMH_init

 A helper function to process the MH_* related initialization.
*********************/
void WtMH_init(WtMHproposal *MHp, 
	     char *MHproposaltype, char *MHproposalpackage, 
	       double *inputs,
	     int fVerbose,
	     WtNetwork *nwp){

  char *fn, *sn;
  int i;
  for (i = 0; MHproposaltype[i] != ' ' && MHproposaltype[i] != 0; i++);
  MHproposaltype[i] = 0;
  /* Extract the required string information from the relevant sources */
  if((fn=(char *)malloc(sizeof(char)*(i+4)))==NULL){
    error("Error in MCMCSample: Can't allocate %d bytes for fn. Memory has not been deallocated, so restart R sometime soon.\n",
	  sizeof(char)*(i+4));
  }
  fn[0]='M';
  fn[1]='H';
  fn[2]='_';
  for(int j=0;j<i;j++)
    fn[j+3]=MHproposaltype[j];
  fn[i+3]='\0';
  /* fn is now the string 'MH_[name]', where [name] is MHproposaltype */
  for (i = 0; MHproposalpackage[i] != ' ' && MHproposalpackage[i] != 0; i++);
  MHproposalpackage[i] = 0;
  if((sn=(char *)malloc(sizeof(char)*(i+1)))==NULL){
    error("Error in ModelInitialize: Can't allocate %d bytes for sn. Memory has not been deallocated, so restart R sometime soon.\n",
	  sizeof(char)*(i+1));
  }
  sn=strncpy(sn,MHproposalpackage,i);
  sn[i]='\0';
  
  /* Search for the MH proposal function pointer */
  MHp->func=(void (*)(WtMHproposal*, WtNetwork*)) R_FindSymbol(fn,sn,NULL);
  if(MHp->func==NULL){
    error("Error in MH_* initialization: could not find function %s in "
	  "namespace for package %s."
	  "Memory has not been deallocated, so restart R sometime soon.\n",fn,sn);
  }

  MHp->inputs=inputs;
  
  MHp->discord=NULL;

  /*Clean up by freeing sn and fn*/
  free((void *)fn);
  free((void *)sn);

  MHp->ntoggles=0;
  (*(MHp->func))(MHp, nwp); /* Call MH proposal function to initialize */
  MHp->toggletail = (Vertex *)malloc(MHp->ntoggles * sizeof(Vertex));
  MHp->togglehead = (Vertex *)malloc(MHp->ntoggles * sizeof(Vertex));
  MHp->toggleweight = (double *)malloc(MHp->ntoggles * sizeof(double));
}

/*********************
 void WtMHProposalInitialize

 A helper function to process the MH_* related initialization.
*********************/
WtMHProposal *WtMHProposalInitialize(
	     char *MHProposaltype, char *MHProposalpackage, 
	       double *inputs,
	     int fVerbose,
	     WtNetwork *nwp){
  WtMHProposal *MHp = Calloc(1, WtMHProposal);

  char *fn, *sn;
  int i;
  for (i = 0; MHProposaltype[i] != ' ' && MHProposaltype[i] != 0; i++);
  MHProposaltype[i] = 0;
  /* Extract the required string information from the relevant sources */
  fn = Calloc(i+4, char);
  fn[0]='M';
  fn[1]='H';
  fn[2]='_';
  for(int j=0;j<i;j++)
    fn[j+3]=MHProposaltype[j];
  fn[i+3]='\0';
  /* fn is now the string 'MH_[name]', where [name] is MHProposaltype */
  for (i = 0; MHProposalpackage[i] != ' ' && MHProposalpackage[i] != 0; i++);
  MHProposalpackage[i] = 0;
  sn = Calloc(i+1, char);
  sn=strncpy(sn,MHProposalpackage,i);
  sn[i]='\0';
  
  /* Search for the MH proposal function pointer */
  MHp->func=(void (*)(WtMHProposal*, WtNetwork*)) R_FindSymbol(fn,sn,NULL);
  if(MHp->func==NULL){
    error("Error in MH_* initialization: could not find function %s in "
	  "namespace for package %s."
	  "Memory has not been deallocated, so restart R sometime soon.\n",fn,sn);
  }

  MHp->inputs=inputs;
  
  MHp->discord=NULL;

  /*Clean up by freeing sn and fn*/
  Free(fn);
  Free(sn);

  MHp->ntoggles=0;
  (*(MHp->func))(MHp, nwp); /* Call MH proposal function to initialize */
  MHp->toggletail = (Vertex *)Calloc(MHp->ntoggles, Vertex);
  MHp->togglehead = (Vertex *)Calloc(MHp->ntoggles, Vertex);
  MHp->toggleweight = (double *)Calloc(MHp->ntoggles, double);

  return MHp;
}

SEXP raw_symbol(SEXP symbol_name) {
  // Find a function pointer for the requested symbol.

  if (!isString(symbol_name) || length(symbol_name)!=1)
    error("Invalid symbol_name.");
  const char *symbol_char = CHAR(STRING_ELT(symbol_name, 0));
  void *symbol = R_FindSymbol(symbol_char, "", NULL);
  if (symbol==NULL)
    error("Could not locate symbol \"%s\".", symbol_char);

  // Copy the function pointer to a raw vector and return it.

  return void_as_raw(symbol);
}

Edge WtEdgeTree2EdgeList(Vertex *tails, Vertex *heads, double *weights, WtNetwork *nwp, Edge nmax){
static Edge (*fun)(Vertex *,Vertex *,double *,WtNetwork *,Edge) = NULL;
if(fun==NULL) fun = (Edge (*)(Vertex *,Vertex *,double *,WtNetwork *,Edge)) R_FindSymbol("WtEdgeTree2EdgeList", "ergm", NULL);
return fun(tails,heads,weights,nwp,nmax);
}

Network * NetworkCopy(Network *src){
static Network * (*fun)(Network *) = NULL;
if(fun==NULL) fun = (Network * (*)(Network *)) R_FindSymbol("NetworkCopy", "ergm", NULL);
return fun(src);
}

void SetEdgeWithTimestamp(Vertex tail, Vertex head, unsigned int weight, Network *nwp){
static void (*fun)(Vertex,Vertex,unsigned int,Network *) = NULL;
if(fun==NULL) fun = (void (*)(Vertex,Vertex,unsigned int,Network *)) R_FindSymbol("SetEdgeWithTimestamp", "ergm", NULL);
fun(tail,head,weight,nwp);
}

void WtNetworkEdgeList(WtNetwork *nwp){
static void (*fun)(WtNetwork *) = NULL;
if(fun==NULL) fun = (void (*)(WtNetwork *)) R_FindSymbol("WtNetworkEdgeList", "ergm", NULL);
fun(nwp);
}

WtModel* WtModelInitialize(char *fnames, char *sonames, double **inputs,int n_terms){
static WtModel* (*fun)(char *,char *,double **,int) = NULL;
if(fun==NULL) fun = (WtModel* (*)(char *,char *,double **,int)) R_FindSymbol("WtModelInitialize", "ergm", NULL);
return fun(fnames,sonames,inputs,n_terms);
}

WtMHProposal * WtMHProposalInitialize(char *MHProposaltype, char *MHProposalpackage,double *inputs,int fVerbose,WtNetwork *nwp){
static WtMHProposal * (*fun)(char *,char *,double *,int,WtNetwork *) = NULL;
if(fun==NULL) fun = (WtMHProposal * (*)(char *,char *,double *,int,WtNetwork *)) R_FindSymbol("WtMHProposalInitialize", "ergm", NULL);
return fun(MHProposaltype,MHProposalpackage,inputs,fVerbose,nwp);
}

/*****************
 int ModelInitialize

 Allocate and initialize the ModelTerm structures, each of which contains
 all necessary information about how to compute one term in the model.
*****************/
Model* ModelInitialize (char *fnames, char *sonames, double **inputsp,
			int n_terms) {
  int i, j, k, l, offset;
  ModelTerm *thisterm;
  char *fn,*sn;
  Model *m;
  double *inputs=*inputsp;
  
  m = (Model *) malloc(sizeof(Model));
  m->n_terms = n_terms;
  m->termarray = (ModelTerm *) malloc(sizeof(ModelTerm) * n_terms);
  m->dstatarray = (double **) malloc(sizeof(double *) * n_terms);
  m->n_stats = 0;
  for (l=0; l < n_terms; l++) {
      thisterm = m->termarray + l;

      /* fnames points to a single character string, consisting of the names
      of the selected options concatenated together and separated by spaces.
      This is passed by the calling R function.  These names are matched with
      their respective C functions that calculate the appropriate statistics. 
      Similarly, sonames points to a character string containing the names
      of the shared object files associated with the respective functions.*/
      for (; *fnames == ' ' || *fnames == 0; fnames++);
      for (i = 0; fnames[i] != ' ' && fnames[i] != 0; i++);
      fnames[i] = 0;
      for (; *sonames == ' ' || *sonames == 0; sonames++);
      for (j = 0; sonames[j] != ' ' && sonames[j] != 0; j++);
      sonames[j] = 0;
      /* Extract the required string information from the relevant sources */
      if((fn=(char *)malloc(sizeof(char)*(i+3)))==NULL){
        error("Error in ModelInitialize: Can't allocate %d bytes for fn. Memory has not been deallocated, so restart R sometime soon.\n",
		sizeof(char)*(i+3));
      }
      fn[0]='d';
      fn[1]='_';
      for(k=0;k<i;k++)
        fn[k+2]=fnames[k];
      fn[i+2]='\0';
      /* fn is now the string 'd_[name]', where [name] is fname */
/*      Rprintf("fn: %s\n",fn); */
      if((sn=(char *)malloc(sizeof(char)*(j+1)))==NULL){
        error("Error in ModelInitialize: Can't allocate %d bytes for sn. Memory has not been deallocated, so restart R sometime soon.\n",
		sizeof(char)*(j+1));
      }
      sn=strncpy(sn,sonames,j);
      sn[j]='\0';
      
      /*  Most important part of the ModelTerm:  A pointer to a
      function that will compute the change in the network statistic of 
      interest for a particular edge toggle.  This function is obtained by
      searching for symbols associated with the object file with prefix
      sn, having the name fn.  Assuming that one is found, we're golden.*/ 
      thisterm->d_func = 
	(void (*)(Edge, Vertex*, Vertex*, ModelTerm*, Network*))
	R_FindSymbol(fn,sn,NULL);
      if(thisterm->d_func==NULL){
        error("Error in ModelInitialize: could not find function %s in "
                "namespace for package %s. Memory has not been deallocated, so restart R sometime soon.\n",fn,sn);
      }      

      /* Optional function to compute the statistic of interest for
	 the network given. It can be more efficient than going one
	 edge at a time. */
      fn[0]='s';
      thisterm->s_func = 
	(void (*)(ModelTerm*, Network*)) R_FindSymbol(fn,sn,NULL);

      /*Clean up by freeing sn and fn*/
      free((void *)fn);
      free((void *)sn);

      /*  Now process the values in model$option[[optionnumber]]$inputs;
          See comments in InitErgm.r for details.    */
      offset = (int) *inputs++;  /* Set offset for attr vector */
/*      Rprintf("offsets: %f %f %f %f %f\n",inputs[0],inputs[1],inputs[2], */
/*		         inputs[3],inputs[4],inputs[5]); */
      thisterm->nstats = (int) *inputs++; /* Set # of statistics returned */
/*      Rprintf("l %d offset %d thisterm %d\n",l,offset,thisterm->nstats); */
      if (thisterm->nstats <= 0)
	{ /* Must return at least one statistic */
	  Rprintf("Error in ModelInitialize:  Option %s cannot return %d \
               statistics.\n", fnames, thisterm->nstats);
	  return NULL;
	}
      /*  Update the running total of statistics */
      m->n_stats += thisterm->nstats; 
      m->dstatarray[l] = (double *) malloc(sizeof(double) * thisterm->nstats);
      thisterm->dstats = m->dstatarray[l];  /* This line is important for
                                               eventually freeing up malloc'ed
					       memory, since thisterm->dstats
					       can be modified but 
					       m->dstatarray[l] cannot be.  */
      thisterm->statcache = (double *) malloc(sizeof(double) * thisterm->nstats);

//.........这里部分代码省略.........

void MCMCSamplePhase12(MHProposal *MH,double *theta, double gain, double *meanstats,int nphase1, int nsubphases, double *networkstatistics,int samplesize, int burnin,int interval, int fVerbose,Network *nwp, Model *m){
static void (*fun)(MHProposal *,double *,double,double *,int,int,double *,int,int,int,int,Network *,Model *) = NULL;
if(fun==NULL) fun = (void (*)(MHProposal *,double *,double,double *,int,int,double *,int,int,int,int,Network *,Model *)) R_FindSymbol("MCMCSamplePhase12", "ergm", NULL);
fun(MH,theta,gain,meanstats,nphase1,nsubphases,networkstatistics,samplesize,burnin,interval,fVerbose,nwp,m);
}

void NetworkDestroy(Network *nwp){
static void (*fun)(Network *) = NULL;
if(fun==NULL) fun = (void (*)(Network *)) R_FindSymbol("NetworkDestroy", "ergm", NULL);
fun(nwp);
}

void MCMCPhase12(int *tails, int *heads, int *dnedges,int *dn, int *dflag, int *bipartite,int *nterms, char **funnames,char **sonames,char **MHProposaltype, char **MHProposalpackage,double *inputs,double *theta0, int *samplesize,double *gain, double *meanstats, int *phase1, int *nsub,double *sample, int *burnin, int *interval,int *newnetworktails,int *newnetworkheads,int *fVerbose,int *attribs, int *maxout, int *maxin, int *minout,int *minin, int *condAllDegExact, int *attriblength,int *maxedges,int *mtails, int *mheads, int *mdnedges){
static void (*fun)(int *,int *,int *,int *,int *,int *,int *,char **,char **,char **,char **,double *,double *,int *,double *,double *,int *,int *,double *,int *,int *,int *,int *,int *,int *,int *,int *,int *,int *,int *,int *,int *,int *,int *,int *) = NULL;
if(fun==NULL) fun = (void (*)(int *,int *,int *,int *,int *,int *,int *,char **,char **,char **,char **,double *,double *,int *,double *,double *,int *,int *,double *,int *,int *,int *,int *,int *,int *,int *,int *,int *,int *,int *,int *,int *,int *,int *,int *)) R_FindSymbol("MCMCPhase12", "ergm", NULL);
fun(tails,heads,dnedges,dn,dflag,bipartite,nterms,funnames,sonames,MHProposaltype,MHProposalpackage,inputs,theta0,samplesize,gain,meanstats,phase1,nsub,sample,burnin,interval,newnetworktails,newnetworkheads,fVerbose,attribs,maxout,maxin,minout,minin,condAllDegExact,attriblength,maxedges,mtails,mheads,mdnedges);
}

void WtChangeStats(unsigned int ntoggles, Vertex *toggletail, Vertex *togglehead, double *toggleweight, WtNetwork *nwp, WtModel *m){
static void (*fun)(unsigned int,Vertex *,Vertex *,double *,WtNetwork *,WtModel *) = NULL;
if(fun==NULL) fun = (void (*)(unsigned int,Vertex *,Vertex *,double *,WtNetwork *,WtModel *)) R_FindSymbol("WtChangeStats", "ergm", NULL);
fun(ntoggles,toggletail,togglehead,toggleweight,nwp,m);
}

void WtModelDestroy(WtModel *m){
static void (*fun)(WtModel *) = NULL;
if(fun==NULL) fun = (void (*)(WtModel *)) R_FindSymbol("WtModelDestroy", "ergm", NULL);
fun(m);
}