# --------- density on right triangle --------------
# Define the locations of the vertices of the right triangle
pts <- matrix(c(0,0, 1,0, 0,1), 3, 2, byrow=TRUE)
npts <- dim(pts)[1]
# Specify the single triangle defined by these vertices.
# The vertices are listed in counter-clockwise order
tri <- matrix(c(1,2,3),1,3)
ntri <- dim(tri)[1]
# Set up the FEM basis object for this region
edg <- NULL
RtTriBasis <- create.FEM.basis(pts, edg, tri, 1)
# List object containing details of nodes.
nodeList <- makenodes(pts,tri,1)
# number of random points required
Nsample <- 100
# generate random points ... define the function
obspts <- randomFEMpts(Nsample, pts, tri)
# plot the random points
plotFEM.mesh(pts, tri)
points(obspts[,1],obspts[,2])
# ----------------- density on a circle ---------------
angle = seq(0,11*pi/6,len=12)
# define the 24 point locations
rad = 2
ctr = c(0,0)
pts = matrix(0,25,2);
pts[ 1:12,1] = rad*cos(angle) + ctr[1]
pts[ 1:12,2] = rad*sin(angle) + ctr[2]
pts[13:24,1] = rad*cos(angle)/2 + ctr[1]
pts[13:24,2] = rad*sin(angle)/2 + ctr[2]
# define the edge indices
edg <- matrix(0,12,2)
for (i in 1:11) {
edg[i,1] <- i
edg[i,2] <- i+1
}
edg[12,] <- c(12,1)
# use delaunay triangulation to define the triangles
tri = delaunayn(pts)
# plot the triangulation of the circle
plotFEM.mesh(pts, tri, xlim=c(-2,2), ylim=c(-2,2))
# Define the true log density function
SinCosIntensFn <- function(x,y) {
return(sin(pi*x/2)*cos(pi*y/2))
}
# locate observation points with sin*cos log density
intDensityVec <- triDensity(pts, tri, SinCosIntensFn)
# number of random points required
Nsample <- 100
# generate random points ... define the function
obspts <- randomFEMpts(Nsample, pts, tri, intDensityVec)
points(obspts[,1], obspts[,2], pch="o")
Run the code above in your browser using DataLab