N = 50; % tamanio de la grilla
V = zeros(N); % matriz rellena con ceros de N x N (tendra los valores de voltajes)
B = -ones(N); % matriz rellena con -1 de N x N (tendra los valores de voltajes en en los bordes)
C = false(N); % matriz binarya de N x N con dos valores 1: donde hay conductor, 0: donde no

% CALCULOS AUXILIARES PARA DEFINIR UN CONDUCTOR CIRCULAR EN EL CENTRO:
x = linspace(-10,10,N); % vector equiespaciado de N elementos que va de -10 a 10
y = linspace(-10,10,N); % vector equiespaciado de N elementos que va de -10 a 10
[X,Y] = meshgrid(x,y); % Construccion de una grilla a partir de vectores (sirve para evaluar funciones de dos variables)
R = sqrt(X.^2 + Y.^2); % R es una funcion que a cada punto del espacio le asigna su distancia al cero
C(R<3) = true;  % Todos los puntos que esten a distancia menor a 3 se ponene a 1 (o true) en la matriz C (representa un conductor circular)

% Definimos el operador que "suaviza"
f = zeros(3);
f(1,2) = 1;
f(2,1) = 1;
f(3,2) = 1;
f(2,3) = 1;
Norm = conv2(ones(N),f,'same'); % Esta matriz sirve para normalizar el suavizado


% Definimos los bordes: 
%B(:,1) = 5;
B(end,:) = 0;
B(1,:) = 5;


for i = 1:10000  % iteramos 10000 veces
V(B>-1) = B(B>-1); % reemplazamos en V los valores de borde
V(C) = mean(V(C)); % promediamos V en toda la region donde hay conductor (sirve solo si hay un unico conductor en la configuracion)

if mod(i,50) == 0  % cada 50 iteraciones graficamos (es solo para tener una idea de como avanza el calculo)
    imagesc(V) %dibujar el estado actual de V
    pause % esperar a tocar una tecla para continuar
end

V = conv2(V,f,'same')./Norm; % ese es el paso cetral (el borroneado)
end