#----------------------------------------------------
#  Campo elettrico di un dipolo
#   dipolo_pm_all.R : mostra contributi delle due cariche
#   (versione che mostra anche il campo fra le cariche)
#
#  GdA 24/10/2021
#----------------------------------------------------

# funzione per calcolare la distanza fra due punti 
# identificati dai vettori r1 e r2
dist <- function(r1,r2) {
   sqrt( (r1[1] - r2[1])^2 + (r1[2] - r2[2])^2 )
}

# funzione per calcolare il modulo di un vettore
mod <- function(v) sqrt( v[1]^2 + v[2]^2 )

k <- 1       # Unità arbitrarie (TUTTE, ma ovviamente consistenti!) 
scala <- 100   # fattore di scala per graficare i vettori forza

a = 6   # distanza fra le cariche
q = 0.1   # modulo della carica

yp = a/2        # posizione carica positiva
ym = -a/2       # posizione carica negativa
rp <- c(0, yp)  # vettori posizione cariche dall'origine
rm <- c(0, ym)

plot(NULL,  xlim=c(-10,10), ylim=c(-10,10), xlab='', ylab='',
     xaxt='n', yaxt='n', asp=1)
abline(h=0, col='gray', lwd=0.5)
abline(v=0, col='gray', lwd=0.5)
abline(0, 1, col='gray', lwd=0.5)
abline(0, -1, col='gray', lwd=0.5)
points(0,yp,pch=19,col='red')
points(0,ym,pch=19,col='black')

cat(sprintf("Clicca per mostrare il campo (nel END per terminare)\n"))
text(9.8,-9.8, 'END', col='gray', cex=2.5)
while(1) {
  p <- locator(1)
  if ( (p$x[1]> 8.1 & p$x[1] < 12) &
       (p$y[1]> -12 & p$y[1] < -9.2) ) break   
  points(p$x[1], p$y[1], pch=19, col='black', cex=0.5)

  r <- c(p$x[1], p$y[1])
  r1 <- r - rp
  r2 <- r - rm
  Ep <- k*q*r1/mod(r1)^3  * scala
  En <- -k*q*r2/mod(r2)^3 * scala

  cat(sprintf("Ep=(%.2e, %.2e); ", Ep[1], Ep[2]))
  cat(sprintf("En=(%.2e, %.2e); ", En[1], En[2]))
  cat(sprintf("E=(%.2e, %.2e) [U.A.]\n", Ep[1]+En[1], Ep[2]+En[2]))
  arrows(r[1], r[2], r[1]+Ep[1], r[2]+Ep[2], length=0.1, col='red')
  arrows(r[1], r[2], r[1]+En[1], r[2]+En[2], length=0.1, col='black')
  arrows(r[1], r[2], r[1]+Ep[1]+En[1], r[2]+Ep[2]+En[2],
         lwd=2, length=0.1, col='blue')
}