Note: If you wish to replicate the R code below, then you will need to copy and paste the following command in R first to make sure you have all the packages and their dependencies installed:
install.packages("install.load", "import", "plot3D", "ggplot2", "rmatio", "reshape2", "directlabels", "gsubfn")
# install the install.load and import packages
## ContourFunctions # contour plot
## use for contour plots -- https://cran.r-project.org/web/packages/ContourFunctions/vignettes/Introduction_to_the_cf_R_package.html
Note: If you wish to replicate the R code below, then you will need to copy and paste the following commands in R first to make sure you have all the packages loaded:
install.load::load_package("plot3D", "ggplot2", "rmatio", "reshape2", "directlabels",
"gsubfn")
# load needed packages using the load_package function from the install.load
# package (it is assumed that you have already installed these packages)
import::from(pracma, meshgrid)
# import meshgrid from the pracma package
The following example, "Force on the Piston of a Pump versus Well Depth and Cylinder Radius", is solved using R (incomplete graphical solution) and GNU Octave. (Brockman 470-473)
The following R code is incomplete as I am still searching for the best way to obtain the same plots in R as with GNU Octave/MATLAB. I posted this question - Q: r - mesh, meshz, surf, and contour plots similar to GNU Octave/MATLAB - online at Stack Overflow, but I have not received an answer yet. Please feel free to post an answer or leave comments. Thank you.
"W = load force"
"r_cylinder = radius of the cylinder, which varies from 15 - 40 mm"
"h = depth of well, which varies from 25 - 50 m"
"rho, density of water, = 1000 kg/m^3"
"W = W(rcylinder, h) = rho * g * (pi * r_cylinder^2 * h)"
Source: Brockman 470 - 473
# R plotmath grDevices - subscripts
# GNU Octave xlabel ylabel subscripts _ underscore
rho <- 1000
g <- 9.81
r_cylinder <- seq(15, 40, by = 5)
h <- seq(25, 50, by = 5)
# mesh plot 'for piston force example data'
list[r_cylinder_grid, h_grid] <- meshgrid(r_cylinder, h)
W_grid <- rho * g * (pi * (r_cylinder_grid/1000)^2 * h_grid)
persp3D(r_cylinder_grid, h_grid, W_grid, xlab = "r_{cylinder}, mm", ylab = "h, m",
zlab = "W, N", theta = 100, phi = -30, expand = 0.5, ticktype = "detailed", ltheta = 120,
bty = "g", facets = FALSE, col = NULL) # Source 1 - 2
# plot_ly(x = r_cylinder, y = h_grid, z = W_grid, type = 'mesh3d')
# meshz plot 'adds 'curtain' around plot to highlight z-values at edges'
# Not available
# surface plot
surf3D(r_cylinder_grid, h_grid, W_grid, xlab = "r_{cylinder}, mm", ylab = "h, m",
zlab = "W, N", theta = 100, phi = -30, expand = 0.5, ticktype = "detailed", ltheta = 120,
bty = "g", col = gg.col(50)) # Source 1
# contour plot
# values obtained from GNU Octave
breaks <- read.mat("C.mat")
breaks <- round(breaks$C, digits = 0)
x <- melt(r_cylinder_grid)
x <- x$value
y <- melt(h_grid)
y <- y$value
z <- melt(W_grid)
z <- z$value
df <- data.frame(x, y, z)
# Source 3 - 4 begins
p <- ggplot(df, aes(x, y, z = z))
p <- p + stat_contour(aes(colour = ..level..), breaks = breaks) + scale_colour_gradient(low = "green",
high = "purple") + theme_bw() + labs(x = bquote(~r[cylinder] ~ ", mm"), y = "h, m",
title = "Contour plot")
direct.label(p, list("far.from.others.borders", "calc.boxes", "enlarge.box", box.color = NA,
fill = "transparent", "draw.rects"))
# Source 3 - 4 ends
# side-view cross-section plots Using base R graphics to plot the matrix plots
# ('side-view cross-section plots')
# 'side view from the 'left'' add r = 15 to r = 40 for h_grid, W_grid
matplot(h_grid, W_grid, type = "l", xlab = "h (m)", ylab = "W (N)")
# 'side view from the 'right'' add h = 25 to h = 50 for r_cylinder_grid, W_grid
matplot(t(r_cylinder_grid), t(W_grid), type = "l", xlab = "r_{cylinder} (mm)", ylab = "W (N)")
rho = 1000;
g = 9.81;
r_cylinder = 15:5:40;
h = 25:5:50;
[r_cylinder_grid, h_grid] = meshgrid(r_cylinder, h);
W_grid = rho*g*(pi*(r_cylinder_grid/1000).^2 .* h_grid);
mesh(r_cylinder_grid, h_grid, W_grid);
xlabel('r_{cylinder}, mm');
ylabel('h, m');
zlabel('W, N');
print("example1.png")
rho = 1000;
g = 9.81;
r_cylinder = 15:5:40;
h = 25:5:50;
[r_cylinder_grid, h_grid] = meshgrid(r_cylinder, h);
W_grid = rho*g*(pi*(r_cylinder_grid/1000).^2 .* h_grid);
meshz(r_cylinder_grid, h_grid, W_grid);
xlabel('r_{cylinder}, mm');
ylabel('h, m');
zlabel('W, N');
print("example2.png")
rho = 1000;
g = 9.81;
r_cylinder = 15:5:40;
h = 25:5:50;
[r_cylinder_grid, h_grid] = meshgrid(r_cylinder, h);
W_grid = rho*g*(pi*(r_cylinder_grid/1000).^2 .* h_grid);
surf(r_cylinder_grid, h_grid, W_grid);
xlabel('r_{cylinder}, mm');
ylabel('h, m');
zlabel('W, N');
print("example3.png")
rho = 1000;
g = 9.81;
r_cylinder = 15:5:40;
h = 25:5:50;
[r_cylinder_grid, h_grid] = meshgrid(r_cylinder, h);
W_grid = rho*g*(pi*(r_cylinder_grid/1000).^2 .* h_grid);
[C, hd1] = contour(r_cylinder_grid, h_grid, W_grid);
clabel(C, hd1);
xlabel('r_{cylinder}, mm');
ylabel('h, m');
% save the object C in the "MATLAB v4 binary data format"
% save("-v4", "C.mat", "C")
print("example4.png")
rho = 1000;
g = 9.81;
r_cylinder = 15:5:40;
h = 25:5:50;
[r_cylinder_grid, h_grid] = meshgrid(r_cylinder, h);
W_grid = rho*g*(pi*(r_cylinder_grid/1000).^2 .* h_grid);
% add r = 15 to r = 40 for h_grid, W_grid
plot(h_grid, W_grid);
xlabel('h, m');
ylabel('W, N');
% add h = 25 to h = 50 for r_cylinder_grid, W_grid
plot(r_cylinder_grid', W_grid')
xlabel('r_{cylinder}, mm');
ylabel('W, N');
print("example5.png")
Source 1
Statistical tools for high-throughput data analysis (STHDA): Impressive package for 3D and 4D graph - R software and data visualization. See http://www.sthda.com/english/wiki/impressive-package-for-3d-and-4d-graph-r-software-and-data-visualization.
Source 2
Fifty ways to draw a volcano using package plot3D by Karline Soetaert. See https://cloud.r-project.org/web/packages/plot3D/vignettes/volcano.pdf.
Source 3
ggplot2 - R: How to Label Specific Contours using direct.label - Stack Overflow answered by a different ben on Jul 21 2014. See https://stackoverflow.com/questions/10675387/r-how-to-label-specific-contours-using-direct-label.
Source 4
R adding legend and directlabels to ggplot2 contour plot - Stack Overflow answered and edited by Sandy Muspratt on Jul 2 2016. See https://stackoverflow.com/questions/38154679/r-adding-legend-and-directlabels-to-ggplot2-contour-plot.
Source 5
r - Superscript and subscript axis labels in ggplot2 - Stack Overflow answered and edited by akrun on Dec 12 2014. See https://stackoverflow.com//questions/27445452/superscript-and-subscript-axis-labels-in-ggplot2.
Source 6
r - RMarkdown Octave Plot - Stack Overflow answered by danlooo on Mar 17, 2022. See https://stackoverflow.com/questions/71512367/rmarkdown-octave-plot
Jay B. Brockman, Introduction to Engineering: Modeling and Problem Solving, Hoboken, New Jersey: John Wiley & Sons, Inc., 2009, pages 470-473.
EcoC2S Home
About EcoC2S
EcoC2S Services
1 Stop Shop
Products
EcoC2S Media
EcoC2S Resources
R Trainings and Resources provided by EcoC2S (Irucka Embry, EIT)
If you would like to contribute to the continued development of Irucka Embry’s R packages and/or Irucka Embry’s R Examples, please feel free to donate via the link below:
Please feel free to review Irucka Embry (iaembry)’s profile on Liberapay.
All R code written by Irucka Embry is distributed under the GPL-3 (or later) license, see the GNU General Public License (GPL) page.
All written content originally created by Irucka Embry is copyrighted under the Creative Commons Attribution-ShareAlike 4.0 International License. All other written content retains the copyright of the original author(s).
This work is licensed under a Creative Commons Attribution-ShareAlike 4.0 International License.