#----------------------------------------------------
#  Campo elettrico di un dipolo
#
#  GdA 23/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 <- 1000   # fattore di scala per graficare i vettori forza

a = 0.1   # 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)
points(0,yp,pch=19,col='red')
points(0,ym,pch=19,col='black')

cat(sprintf("Clicca per generare vettori di campo (sul dipolo per terminare)\n"))
while(1) {
  p <- locator(1)
  if ( mod(c(p$x[1],p$y[1])) < a*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("E=(%.2e, %.2e) [U.A.]\n", Ep[1]+En[1], Ep[2]+En[2]))
  arrows(r[1], r[2], r[1]+Ep[1]+En[1], r[2]+Ep[2]+En[2],
         lwd=2, length=0.1, col='blue')
}