rm(list=ls())
#
# Funktion zur Berechnung des Schnittpunktes zwischen Linie und Kreis
#
int.sect <- function(rad,cntr,c,d) {
  A <- d^2+1
  B <- 2*c*d - 2*cntr
  C <- rad^2 - cntr^2 -c^2
  D <- C/A + B^2/(4*A^2)
  x1 <- -B/(2*A) + sqrt(D) 
  x2 <- -B/(2*A) - sqrt(D)
  y1 <- c + d*x1
  y2 <- c + d*x2
  int.sect <- c(x1,x2,y1,y2,D)
} 
#
# --- Der Mohr'sche Spannungskreis ---
#
sigma3 <-  50  # kN/m^2    # 1. Hauptspannung
sigma1 <- 150  # kN/m^2   # 3. Hauptspannung
#
# Der Kreis
#
narc <- 150
Arc <- seq(0,pi/2,pi/(2*narc))
fun_sigma <- function(sigma1,sigma3,Arc){
                0.5*((sigma1+sigma3)+(sigma1-sigma3)*cos(2*Arc))
             }
fun_tau   <- function(sigma1,sigma3,Arc){
                0.5*(sigma1-sigma3)*sin(2*Arc)
             }
tau   <- fun_tau(sigma1,sigma3,Arc)
sigma <- fun_sigma(sigma1,sigma3,Arc)
plot(sigma,tau,type="l",xlim=c(0,1.1*max(sigma)),asp=1,lwd=2,
       ylim=c(-5,1.2*max(tau)),pch=20,cex=0.8,col="darkgreen",
       main="Mohr-Coulomb'sche\nGrenzgerade",
       xlab=expression(paste("Normalspannung  ",sigma)),
       ylab=expression(paste("Scherspannung   ",tau)))
text(sigma[1],-5,expression(sigma[1]))
text(sigma[narc],-5,expression(sigma[3]))
lines(c(-sigma1,2*sigma1),c(0,0))
#
# --- Die MC-Grenzgerade ---
#
coh <- 25
arc.fric <- 14.9808
# Test auf Schnittpunkte:
tan.fric <- tan(pi*arc.fric/180) 
rad <- 0.5*abs(sigma1-sigma3)
cntr <- 0.5*(sigma1 + sigma3)
int.point <- int.sect(rad,cntr,coh,tan.fric)
# Zeichnen der Ergebnisse:
if(int.point[5]>=0) 
{
  lcl <- "red"
  lines(c(0,max(sigma1,sigma3)),c(coh,coh+tan.fric*max(sigma1,sigma3)),
        col=lcl,lwd=2)
  points(int.point[1:2],int.point[3:4])
} else {        
  lcl <- "blue"
  lines(c(0,max(sigma1,sigma3)),c(coh,coh+tan.fric*max(sigma1,sigma3)),
        col=lcl,lwd=2)
}