rm(list=ls()) 
# Load required packages
library(readr)
library(ggplot2)
library(tidyr)
library(plotly)
library(processx)
library(reshape2)
#font
t <- list(
  family = "Arial",
  size = 14,
  color = 'black')

relative density vs C.

# Read the CSV file
gyroid_u1.2<-read.csv(file='gyroid_geometry_cong_u1.2.csv', header = TRUE, dec='.', sep=',')

# rela_density ~ C: rela_density = a*C + b
lm_model <- lm(rela_density ~ C, data = gyroid_u1.2)
# Extract the coefficients
a <- coef(lm_model)[2]  # coefficient for C
b <- coef(lm_model)[1]  # intercept term
# Print the coefficients
cat("Coefficient a:", a, "\n")
## Coefficient a: -0.3331213
cat("Coefficient b:", b, "\n")
## Coefficient b: 0.5000124
# Get the model summary
summary(lm_model)
## 
## Call:
## lm(formula = rela_density ~ C, data = gyroid_u1.2)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -0.0119839 -0.0039080  0.0001466  0.0038832  0.0119592 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.5000124  0.0009114   548.6   <2e-16 ***
## C           -0.3331213  0.0010499  -317.3   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.004984 on 28 degrees of freedom
## Multiple R-squared:  0.9997, Adjusted R-squared:  0.9997 
## F-statistic: 1.007e+05 on 1 and 28 DF,  p-value: < 2.2e-16
# Calculate the correlation coefficient
correlation <- cor(gyroid_u1.2$rela_density, gyroid_u1.2$C)
cat("Correlation coefficient:", correlation, "\n")
## Correlation coefficient: -0.999861
# Create scatter plot
scatter_plot_1 <- ggplot(gyroid_u1.2, aes(x = C, y = rela_density)) +
  geom_point() +
  labs(x = "C", y = "Relative density [%]") +
  scale_x_continuous(breaks = seq(-1.5, 1.5, 0.6)) +
  scale_y_continuous(breaks = seq(0, 1, 0.2), labels = paste0(seq(0, 100, 20), "%")) +
  
    theme(
      axis.text = element_text(family = "Arial",size = 10,color = "black"),
      axis.title = element_text(family = "Arial",size = 10),
      axis.ticks = element_line(color = "black", size = 0.5),
      axis.ticks.length = unit(-0.1, "mm"),
      panel.background = element_rect(fill = "white", color = "black"),
      panel.border = element_rect(color = "black",size = 1, fill = NA),
      #panel.grid.major.y = element_line(color = "gray80"),  # Set the major grid lines color
      #panel.grid.major.x = element_line(color = "gray80"),
    )+
  geom_abline(intercept = b, slope = a, linetype = "dashed", color = "black")
scatter_plot_1

feature size in Gyroid structures

pore_channel, pore_core, strut diameter shows “Hyperbolic paraboloid” with unit cell size and C (relative density). analytical result as red surfaces

# Read the CSV file
gyroid<-read.csv(file='gyroid_geometry_cong_c-1to1.csv', header = TRUE, dec='.', sep=',')
gyroid$unit<-as.numeric(gyroid$unit)
gyroid_dia<-read.csv(file='gyroid_geometry_cong_c-1to1_dia.csv', header = TRUE, dec='.', sep=',')
gyroid_dia$unit<-as.numeric(gyroid_dia$unit)

3d display for pore core (unit and C)

pore_fig<-plot_ly(gyroid,x = ~unit, y = ~C, z = ~pore_core,
  marker = list(color = ~pore_core, colorscale = 'Greens', cmin = -0.2, cmax = 3,showscale = TRUE)) %>%add_markers()



pore_fig <- pore_fig %>% layout(scene = list(xaxis = list(title = 'Unit cell size [mm]',
                                                     tickvals =c(1, 2, 3, 4),
                                                ticktext = paste0(c(1, 2, 3, 4))),
                                   yaxis = list(title = 'C',
                                                     tickvals =c(-1.0,-0.8,-0.5,-0.3,0.0,0.3,0.5,0.8,1.0),
                                                ticktext = paste0(c(-1.0,"",-0.5,"",0.0,"",0.5,"",1.0))),
                                   zaxis = list(title = 'Pore core size [mm]',
                                                range = c(-0.2, 3.0)),
                                   #camera perspective
                                   camera = list(eye = list(x = -1.25, y = -1.25, z = 0.1))),
                                font=t,
                      annotations = list(
                        x = 1.13,
                        y = 1.05,
                        z = 1,
                        font=t,
                        text = 'Pore core size',
                        xref = 'paper',
                        yref = 'paper',
                        showarrow = FALSE
                        ))
# Create the plane based on the equation
dummy_unit <- seq(1, 4, length.out = 9)
dummy_c <- seq(-0.9, 1, length.out = 20)
#an analytical surface for pore
ana_pore_core_data<-expand.grid(unit = dummy_unit,C = dummy_c,KEEP.OUT.ATTRS = F)
ana_pore_core_data$pore_core <-ana_pore_core_data$unit+ana_pore_core_data$unit*(atan((1-sqrt(2-ana_pore_core_data$C*ana_pore_core_data$C))/(ana_pore_core_data$C+1)))/pi-ana_pore_core_data$unit*(atan((1+sqrt(2-ana_pore_core_data$C*ana_pore_core_data$C))/(ana_pore_core_data$C+1)))/pi
ana_pore_core_data <- acast(ana_pore_core_data, C~unit , value.var = "pore_core") #y ~ x

# add analytical surface into 3d plot
pore_fig <- add_trace(p=pore_fig, type = "surface",
                      x = dummy_unit,
                      y = dummy_c,
                      z = ana_pore_core_data,
                      type = "surface",
                      showscale = FALSE,
                      colorscale = list(c(0, 1), c("#f0202080", "#f0202080"))
                      )
#
pore_fig

3d display for pore channel (unit and C)

pore_fig_2<-plot_ly(gyroid, x = ~unit, y = ~C, z = ~pore_channel,
                  marker = list(color = ~pore_channel,  colorscale = 'Greens', cmin = -0.2, cmax = 3,showscale = TRUE)) %>% add_markers()

pore_fig_2 <- pore_fig_2 %>% layout(scene = list(xaxis = list(title = 'Unit cell size [mm]',
                                                     tickvals =c(1, 2, 3, 4),
                                                ticktext = paste0(c(1, 2, 3, 4))),
                                   yaxis = list(title = 'C',
                                                     tickvals =c(-1.0,-0.8,-0.5,-0.3,0.0,0.3,0.5,0.8,1.0),
                                                ticktext = paste0(c(-1.0,"",-0.5,"",0.0,"",0.5,"",1.0))),
                                   zaxis = list(title = 'Pore channel size [mm]',
                                                range = c(-0.2, 3.0)),
                                   #camera perspective
                                   camera = list(eye = list(x = -1.25, y = -1.25, z = 0.1))),
                                font=t,
                      annotations = list(
                        x = 1.13,
                        y = 1.05,
                        z = 1,
                        font=t,
                        text = 'Pore channel size',
                        xref = 'paper',
                        yref = 'paper',
                        showarrow = FALSE
                        ))
# Create the plane based on the equation
#axis
range_dummy_unit <- range(gyroid$unit)
range_dummy_c <- range(gyroid$C)
dummy_unit <- seq(1, 4, length.out = 9)
dummy_c <- seq(-0.9, 1, length.out = 20)

#an analytical surface for pore
ana_pore_channel_data<-expand.grid(unit = dummy_unit,C = dummy_c,KEEP.OUT.ATTRS = F)
ana_pore_channel_data$pore_channel <-ana_pore_channel_data$unit-2*ana_pore_channel_data$unit*(atan((1+sqrt(2-ana_pore_channel_data$C*ana_pore_channel_data$C))/(ana_pore_channel_data$C+1)))/pi
ana_pore_channel_data <- acast(ana_pore_channel_data, C~unit , value.var = "pore_channel") #y ~ x


pore_fig_2 <- add_trace(p=pore_fig_2, type = "surface",
                      x = dummy_unit,
                      y = dummy_c,
                      z = ana_pore_channel_data,
                      type = "surface",
                      showscale = FALSE,
                      colorscale = list(c(0, 1), c("#f0202080", "#f0202080"))
                      )
#
pore_fig_2

3d display for strut diameter (unit and C)

pore_fig_3<-plot_ly(gyroid_dia, x = ~unit,  y = ~C,  z = ~strut_dia,
                  marker = list(color = ~strut_dia, colorscale = 'Greens', cmin = -0.2, cmax = 3, showscale = TRUE)) %>% add_markers()

pore_fig_3 <- pore_fig_3 %>% layout(scene = list(xaxis = list(title = 'Unit cell size [mm]',
                                                     tickvals =c(1, 2, 3, 4),
                                                ticktext = paste0(c(1, 2, 3, 4))),
                                   yaxis = list(title = 'C',
                                                     tickvals =c(-1.0,-0.8,-0.5,-0.3,0.0,0.3,0.5,0.8,1.0),
                                                ticktext = paste0(c(-1.0,"",-0.5,"",0.0,"",0.5,"",1.0))),
                                   zaxis = list(title = 'Strut diameter [mm]',
                                                range = c(-0.2, 3.0)),
                                   #camera perspective
                                   camera = list(eye = list(x = -1.25, y = -1.25, z = 0.1))),
                                font=t,
                      annotations = list(
                        x = 1.13,
                        y = 1.05,
                        z = 1,
                        font=t,
                        text = 'Strut diameter',
                        xref = 'paper',
                        yref = 'paper',
                        showarrow = FALSE
                        ))
# Create the plane based on the equation
#axis
range_dummy_unit <- range(gyroid_dia$unit)
range_dummy_c <- range(gyroid_dia$C)
dummy_unit <- seq(1, 4, length.out = 9)
dummy_c <- seq(-0.9, 1, length.out = 20)
#an analytical surface for pore
ana_strut_dia_data<-expand.grid(unit = dummy_unit,C = dummy_c,KEEP.OUT.ATTRS = F)
ana_strut_dia_data$strut_dia<-sqrt(2)*0.5*ana_strut_dia_data$unit-ana_strut_dia_data$unit+2*ana_strut_dia_data$unit*(atan((1+sqrt(2-ana_strut_dia_data$C*ana_strut_dia_data$C))/(ana_strut_dia_data$C+1)))/pi
ana_strut_dia_data <- acast(ana_strut_dia_data, C~unit , value.var = "strut_dia") #y ~ x

# add plane into 3d plot

pore_fig_3 <- add_trace(p=pore_fig_3, type = "surface",
                      x = dummy_unit, 
                      y = dummy_c, 
                      z = ana_strut_dia_data,
                      type = "surface",
                      showscale = FALSE,
                      colorscale = list(c(0, 1), c("#f0202080", "#f0202080"))
                      )
#
pore_fig_3