Problem:
 Score : 700 points    Problem Statement Nuske has a grid with N rows and M columns of squares. The rows are numbered 1 through N from top to bottom, and the columns are numbered 1 through M from left to right.
Each square in the grid is painted in either blue or white. If S_{i,j} is 1 , the square at the i -th row and j -th column is blue; if S_{i,j} is 0 , the square is white.
For every pair of two blue square a and b , there is at most one path that starts from a , repeatedly proceeds to an adjacent (side by side) blue square and finally reaches b , without traversing the same square more than once.  Phantom Thnook, Nuske's eternal rival, gives Q queries to Nuske. The i -th query consists of four integers x_{i,1} , y_{i,1} , x_{i,2} and y_{i,2} and asks him the following: when the rectangular region of the grid bounded by (and including) the x_{i,1} -th row, x_{i,2} -th row, y_{i,1} -th column and y_{i,2} -th column is cut out, how many connected components consisting of blue squares there are in the region?  Process all the queries.      Constraints  1 ≤ N,M ≤ 2000  1 ≤ Q ≤ 200000  S_{i,j} is either 0 or 1 .  S_{i,j} satisfies the condition explained in the statement.  1 ≤ x_{i,1} ≤ x_{i,2} ≤ N(1 ≤ i ≤ Q)  1 ≤ y_{i,1} ≤ y_{i,2} ≤ M(1 ≤ i ≤ Q)         Input The input is given from Standard Input in the following format:  N  M  Q  S_{1,1} .. S_{1,M} : S_{N,1} .. S_{N,M}  x_{1,1}  y_{i,1}  x_{i,2}  y_{i,2} : x_{Q,1}  y_{Q,1}  x_{Q,2}  y_{Q,2}       Output For each query, print the number of the connected components consisting of blue squares in the region.        Sample Input 1 3 4 4
1101
0110
1101
1 1 3 4
1 1 3 1
2 2 3 4
1 2 2 4      Sample Output 1 3
2
2
2   In the first query, the whole grid is specified. There are three components consisting of blue squares, and thus 3 should be printed.  In the second query, the region within the red frame is specified. There are two components consisting of blue squares, and thus 2 should be printed.
Note that squares that belong to the same component in the original grid may belong to different components.       Sample Input 2 5 5 6
11010
01110
10101
11101
01010
1 1 5 5
1 2 4 5
2 3 3 4
3 3 3 3
3 1 3 5
1 1 3 4      Sample Output 2 3
2
1
1
3
2

Solution (Python):
def main():
    import sys
    import numpy as np
    input = sys.stdin.readline
    H, W, Q = map(int, input().split())
    S = [list(map(lambda x: x=="1", input())) for _ in range(H)] + [[False]*(W+1)]
    S = np.array(S, dtype=np.int16)
    M = np.zeros((H + 1, W + 1), dtype=np.int32)

    UD = S[1:, :] * S[:-1, :]
    M[1:, :] += UD
    M[:-1, :] += UD

    LR = S[:, 1:] * S[:, :-1]
    M[:, 1:] += LR
    M[:, :-1] += LR

    M = S * 2 - M
    M = M.cumsum(axis=0).cumsum(axis=1)
    UD = UD.cumsum(axis=1, dtype=np.int16)
    LR = LR.cumsum(axis=0, dtype=np.int16)
    M[:, -1] = 0
    M[-1, :] = 0

    del S

    UD = np.pad(UD, ((0, 1), (0, 1)), mode="constant")
    LR = np.pad(LR, ((0, 1), (0, 1)), mode="constant")
    x1, y1, x2, y2 = np.array(list(map(int, sys.stdin.read().split())), dtype=np.int16).reshape(-1, 4).T
    Ans = M[x2-1, y2-1] + M[x1-2, y1-2] - M[x2-1, y1-2] - M[x1-2, y2-1]
    Ans += LR[x2-1, y1-2] - LR[x1-2, y1-2]
    Ans += LR[x2-1, y2-1] - LR[x1-2, y2-1]
    Ans += UD[x1-2, y2-1] - UD[x1-2, y1-2]
    Ans += UD[x2-1, y2-1] - UD[x2-1, y1-2]
    Ans //= 2
    print("\n".join(str(int(a)) for a in Ans))

main()

Solution (Python):
def main():
    import sys
    input = sys.stdin.readline
    H, W, Q = map(int, input().split())
    S = [list(map(lambda x: x=="1", input())) for _ in range(H)]
    

    M, L, R, U, D = [], [], [], [], []
    for x, s in enumerate(S):
        M_, L_, R_, U_, D_ = [], [], [], [], []
        for y, c in enumerate(s):
            m = l = r = u = d = 0
            if c:
                if x >= 1:
                    m += S[x-1][y]
                    u -= S[x-1][y]
                if x < H-1:
                    m += S[x+1][y]
                    d -= S[x+1][y]
                if y >= 1:
                    m += s[y-1]
                    l -= s[y-1]
                if y < W-1:
                    m += s[y+1]
                    r -= s[y+1]
                m = 2 - m
                
                
                
                
            M_.append(m)
            L_.append(l)
            R_.append(r)
            U_.append(u)
            D_.append(d)
        M.append(M_)
        L.append(L_)
        R.append(R_)
        U.append(U_)
        D.append(D_)
    M.append([0]*(W+1))
    L.append([0]*(W+1))
    R.append([0]*(W+1))
    U.append([0]*(W+1))
    D.append([0]*(W+1))

    
    
    for x in range(H):
        m = 0
        for y in range(W):
            m += M[x][y]
            M[x][y] = m + M[x-1][y]
            L[x][y] += L[x-1][y]
            R[x][y] += R[x-1][y]
            U[x][y] += U[x][y-1]
            D[x][y] += D[x][y-1]

    Ans = []
    
    for _ in range(Q):
        x1, y1, x2, y2 = map(int, input().split())
        ans = 0
        ans += M[x2-1][y2-1] + M[x1-2][y1-2] - M[x2-1][y1-2] - M[x1-2][y2-1]
        ans -= L[x2-1][y1-1] - L[x1-2][y1-1]
        ans -= R[x2-1][y2-1] - R[x1-2][y2-1]
        ans -= U[x1-1][y2-1] - U[x1-1][y1-2]
        ans -= D[x2-1][y2-1] - D[x2-1][y1-2]
        assert ans % 2 == 0, ans
        
        print(ans//2)
    

main()

Solution (Python):
import sys
def input():
	return sys.stdin.readline()[:-1]

n, m, q = map(int, input().split())
s = ["0"*(m+1)] + ["0" + input() for _ in range(n)]

block = [[0 for _ in range(m+1)] for _ in range(n+1)]
for i in range(1, n+1):
	for j in range(1, m+1):
		block[i][j] = block[i][j-1] + block[i-1][j] - block[i-1][j-1] + int(s[i][j])

def block_num(x1, y1, x2, y2):
	return block[x2][y2] - block[x1-1][y2] - block[x2][y1-1] + block[x1-1][y1-1]

edge_h = [[0 for _ in range(m)] for _ in range(n+1)]
edge_v = [[0 for _ in range(m+1)] for _ in range(n)]

for i in range(1, n+1):
	for j in range(1, m):
		edge_h[i][j] = edge_h[i][j-1] + edge_h[i-1][j] - edge_h[i-1][j-1] + (int(s[i][j]) + int(s[i][j+1]) == 2)

for i in range(1, n):
	for j in range(1, m+1):
		edge_v[i][j] = edge_v[i][j-1] + edge_v[i-1][j] - edge_v[i-1][j-1] + (int(s[i][j]) + int(s[i+1][j]) == 2)

def edge_num(x1, y1, x2, y2):
	h = edge_h[x2][y2-1] - edge_h[x1-1][y2-1] - edge_h[x2][y1-1] + edge_h[x1-1][y1-1]
	v = edge_v[x2-1][y2] - edge_v[x1-1][y2] - edge_v[x2-1][y1-1] + edge_v[x1-1][y1-1]
	return h+v

for _ in range(q):
	x1, y1, x2, y2 = map(int, input().split())
	b = block_num(x1, y1, x2, y2)
	e = edge_num(x1, y1, x2, y2)
	print(b-e)

Solution (Python):
import sys
input = sys.stdin.readline

import numpy as np

N,M,Q = map(int,sys.stdin.readline().split())
rows = sys.stdin.readlines()
S = np.zeros((N+2,M+2),dtype='U1')
S[1:-1,1:-1] = np.array([list(row) for row in rows[:N]],dtype='U1')[:,:M]
query = np.array([row.split() for row in rows[N:]],dtype=np.int32)

count_one = (S=='1').cumsum(axis=0).cumsum(axis=1)


count_edge_x = ((S[:-1]=='1')&(S[1:]=='1')).cumsum(axis=0).cumsum(axis=1)
count_edge_y = ((S[:,:-1]=='1')&(S[:,1:]=='1')).cumsum(axis=0).cumsum(axis=1)

x1,y1,x2,y2 = query.T

answer = count_one[x2,y2] - count_one[x2,y1-1] - count_one[x1-1,y2] + count_one[x1-1,y1-1]

answer -= count_edge_x[x2-1,y2] - count_edge_x[x2-1,y1-1] - count_edge_x[x1-1,y2] + count_edge_x[x1-1,y1-1]
answer -= count_edge_y[x2,y2-1] - count_edge_y[x1-1,y2-1] - count_edge_y[x2,y1-1] + count_edge_y[x1-1,y1-1]

print('\n'.join(answer.astype(str)))

Solution (Python):
N,M,Q=map(int,input().split())
G=[input() for i in range(N)]
node=[[0]*(M+1)]
edge_x=[[0]*M]
edge_y=[[0]*(M+1)]
for i in range(1,N+1):
  n,e_x,e_y=[0],[0],[0]
  for j in range(1,M+1):
    if G[i-1][j-1]=="1":
      n.append(n[-1]+node[i-1][j]-node[i-1][j-1]+1)
    else:
      n.append(n[-1]+node[i-1][j]-node[i-1][j-1])
    if j<M:
      if G[i-1][j-1]=="1" and G[i-1][j]=="1":
        e_x.append(e_x[-1]+edge_x[i-1][j]-edge_x[i-1][j-1]+1)
      else:
        e_x.append(e_x[-1]+edge_x[i-1][j]-edge_x[i-1][j-1])
    if i<N:
      if G[i-1][j-1]=="1" and G[i][j-1]=="1":
        e_y.append(e_y[-1]+edge_y[i-1][j]-edge_y[i-1][j-1]+1)
      else:
        e_y.append(e_y[-1]+edge_y[i-1][j]-edge_y[i-1][j-1])
  
  node.append(n)
  edge_x.append(e_x)
  if i<N:
    edge_y.append(e_y)

for q in range(Q):
  x1,y1,x2,y2=map(int,input().split())
  n=node[x2][y2]-node[x1-1][y2]-node[x2][y1-1]+node[x1-1][y1-1]
  e_y=edge_y[x2-1][y2]-edge_y[x1-1][y2]-edge_y[x2-1][y1-1]+edge_y[x1-1][y1-1]
  e_x=edge_x[x2][y2-1]-edge_x[x1-1][y2-1]-edge_x[x2][y1-1]+edge_x[x1-1][y1-1]
  e=e_x+e_y
  print(n-e)

Solution (Python):
N,M,Q = map(int,input().split())
src = [input() for i in range(N)]
qs = [tuple(map(int,input().split())) for i in range(Q)]

blue_cum = [[0 for j in range(M+1)] for i in range(N+1)]
for i in range(N):
    for j in range(M):
        blue_cum[i+1][j+1] += blue_cum[i+1][j] + int(src[i][j])
for i in range(N):
    for j in range(M+1):
        blue_cum[i+1][j] += blue_cum[i][j]

h_edge = [[0 for j in range(M)] for i in range(N+1)]
for i in range(N):
    for j in range(M-1):
        h_edge[i+1][j+1] += h_edge[i+1][j] + int(src[i][j:j+2]=='11')
for i in range(N):
    for j in range(M):
        h_edge[i+1][j] += h_edge[i][j]

v_edge = [[0 for j in range(M+1)] for i in range(N)]
for i in range(N-1):
    for j in range(M):
        v_edge[i+1][j+1] += v_edge[i+1][j] + int(src[i][j]=='1' and src[i+1][j]=='1')
for i in range(N-1):
    for j in range(M+1):
        v_edge[i+1][j] += v_edge[i][j]

ans = []
for y1,x1,y2,x2 in qs:
    y1 -= 1
    x1 -= 1
    blue = blue_cum[y1][x1] - blue_cum[y1][x2] - blue_cum[y2][x1] + blue_cum[y2][x2]
    edge = h_edge[y1][x1] - h_edge[y1][x2-1] - h_edge[y2][x1] + h_edge[y2][x2-1]
    edge+= v_edge[y1][x1] - v_edge[y1][x2] - v_edge[y2-1][x1] + v_edge[y2-1][x2]
    ans.append(blue - edge)
print(*ans, sep='\n')

Solution (Python):
f=lambda s,h,w,a,q:s[h+1].__setitem__(w+1,s[h+1][w]+s[h][w+1]-s[h][w]+a[h][w]*q);g=lambda s,y1,x1,y2,x2:s[y2][x2]-s[y2][x1-1]-s[y1-1][x2]+s[y1-1][x1-1];H,W,N=map(int,raw_input().split());a=[map(int,raw_input())for _ in[0]*H];s,s1,s2=[[[0]*(W+1)for _ in[0]*-~H]for _ in[0]*3]
for h in range(H):
	for w in range(W):f(s,h,w,a,1);f(s1,h,w,a,w+1<W and a[h][w+1]);f(s2,h,w,a,h+1<H and a[h+1][w])
for _ in[0]*N:y1,x1,y2,x2=map(int,raw_input().split());print g(s,y1,x1,y2,x2)-g(s1,y1,x1,y2,x2-1)-g(s2,y1,x1,y2-1,x2)

Solution (Python):
import math,string,itertools,fractions,heapq,collections,re,array,bisect,sys,random,time,copy,functools

sys.setrecursionlimit(10**7)
inf = 10**20
mod = 10**9 + 7

def LI(): return [int(x) for x in sys.stdin.readline().split()]
def LI_(): return [int(x)-1 for x in sys.stdin.readline().split()]
def LF(): return [float(x) for x in sys.stdin.readline().split()]
def LS(): return sys.stdin.readline().split()
def I(): return int(sys.stdin.readline())
def F(): return float(sys.stdin.readline())
def S(): return input()

class Ruiwa():
    def __init__(self, a):
        self.H = h = len(a)
        self.W = w = len(a[0])
        self.R = r = a
        for i in range(h):
            for j in range(1,w):
                r[i][j] += r[i][j-1]

        for i in range(1,h):
            for j in range(w):
                r[i][j] += r[i-1][j]

    def search(self, x1, y1, x2, y2):
        if x1 > x2 or y1 > y2:
            return 0

        r = self.R
        rr = r[y2][x2]
        if x1 > 0 and y1 > 0:
            return rr - r[y1-1][x2] - r[y2][x1-1] + r[y1-1][x1-1]
        if x1 > 0:
            rr -= r[y2][x1-1]
        if y1 > 0:
            rr -= r[y1-1][x2]

        return rr


def main():
    n,m,q = LI()
    a = [[1 if c == '1' else 0 for c in S()] for _ in range(n)]
    b = [[0] * m for _ in range(n)]
    c = [[0] * m for _ in range(n)]
    for i in range(n):
        for j in range(m):
            if a[i][j] == 0:
                continue
            if i > 0 and a[i-1][j]:
                b[i][j] += 1
            if j > 0 and a[i][j-1]:
                c[i][j] += 1
    A = Ruiwa(a)
    B = Ruiwa(b)
    C = Ruiwa(c)
    rr = []
    for i in range(q):
        y1,x1,y2,x2 = LI_()
        r = A.search(x1,y1,x2,y2) - B.search(x1,y1+1,x2,y2) - C.search(x1+1,y1,x2,y2)
        rr.append(r)

    return '\n'.join(map(str,rr))


print(main())

Solution (Python):
n,m,q=map(int,input().split())
s,v,el,eu=([[0]*(m+2) for _ in range(n+2)] for _ in range(4))
for i in range(1,n+1):
    s[i][1:m+1]=map(int,input())
    for j in range(1,m+1):
        v[i][j]=v[i][j-1]+s[i][j]
        el[i][j]=el[i][j-1]+s[i][j]*s[i][j-1]
        eu[i][j]=eu[i][j-1]+s[i][j]*s[i-1][j]
for i in range(2,n+1):
    for j in range(1,m+1):
        v[i][j]+=v[i-1][j]
        el[i][j]+=el[i-1][j]
        eu[i][j]+=eu[i-1][j]
for _ in range(q):
    x1,y1,x2,y2=map(int,input().split())
    print(v[x2][y2]-v[x1-1][y2]-v[x2][y1-1]+v[x1-1][y1-1]
            -(el[x2][y2]-el[x1-1][y2]-el[x2][y1]+el[x1-1][y1])
            -(eu[x2][y2]-eu[x1][y2]-eu[x2][y1-1]+eu[x1][y1-1]))

Solution (Python):
A,B=raw_input,range;n,m,q=map(int,A().split());s=[map(int,A())for _ in[0]*n]
def f(a):
 z=[[x*y for(x,y)in zip(t,u)]+[0,0]for(t,u)in zip(s,a)]+[[0]*m,[0]*m]
 def g(p):
	for i in B(n):
	 for j in B(m):z[i][j]+=z[i-p][j-1+p]
 g(1);g(0);return z
a,b,c=f([[1]*m]*n),f(t[1:]for t in s),f(s[1:])
for _ in[0]*q:x,y,z,w=map(int,A().split());g=lambda d,p,q:d[x-2][y-2]+d[z-p][w-q]-d[x-2][w-q]-d[z-p][y-2];print g(a,1,1)-g(b,1,2)-g(c,2,1)

Solution (Python):
import sys
from itertools import accumulate

def solve():
    n, m, q = map(int, input().split())
    s = [[int(j) for j in input()] for i in range(n)]

    
    sxy = [list(accumulate(s[i])) for i in range(n)]

    for j in range(m):
        for i in range(n - 1):
            sxy[i + 1][j] += sxy[i][j]

    
    exy = [[0]*(2*m) for i in range(2*n)]

    for i in range(2*n - 1):
        if i % 2 == 0:
            for j in range(2*m - 1):
                exy[i][j + 1] += exy[i][j]
                if j % 2 == 0 and j < 2*m - 2:
                    exy[i][j + 1] += (s[i//2][j//2] == s[i//2][j//2 + 1] == 1)
        else:
            for j in range(2*m - 1):
                if j % 2 == 0:
                    exy[i][j] += (s[i//2][j//2] == s[i//2 + 1][j//2] == 1)

                exy[i][j + 1] += exy[i][j]

    for j in range(2*m):
        for i in range(2*n - 1):
            exy[i + 1][j] += exy[i][j]

    '''
    for i in range(2*n):
        print(exy[i], file=sys.stderr)
    '''

    for qi in range(q):
        x1, y1, x2, y2 = map(lambda x: x - 1, map(int, sys.stdin.readline().split()))

        nv = sxy[x2][y2]

        if x1 - 1 >= 0:
            nv -= sxy[x1 - 1][y2]
        if y1 - 1 >= 0:
            nv -= sxy[x2][y1 - 1]
        if x1 - 1 >= 0 and y1 - 1 >= 0:
            nv += sxy[x1 - 1][y1 - 1]

        ne = exy[2*x2][2*y2]

        if 2*x1 - 1 >= 0:
            ne -= exy[2*x1 - 1][2*y2]
        if 2*y1 - 1 >= 0:
            ne -= exy[2*x2][2*y1 - 1]
        if 2*x1 - 1 >= 0 and 2*y1 - 1 >= 0:
            ne += exy[2*x1 - 1][2*y1 - 1]

        ans = nv - ne

        print(ans)

if __name__ == '__main__':
    solve()

Solution (Python):
f=lambda s,h,w,a,q:s[h+1][w]+s[h][w+1]-s[h][w]+a[h][w]*q;g=lambda s,y1,x1,y2,x2:s[y2][x2]-s[y2][x1-1]-s[y1-1][x2]+s[y1-1][x1-1];H,W,N=map(int,raw_input().split());a=[map(int,raw_input())for _ in[0]*H];s,s1,s2=[[[0]*(W+1)for _ in[0]*-~H]for _ in[0]*3]
for h in range(H):
	for w in range(W):s[h+1][w+1]=f(s,h,w,a,1);s1[h+1][w+1]=f(s1,h,w,a,w+1<W and a[h][w+1]);s2[h+1][w+1]=f(s2,h,w,a,h+1<H and a[h+1][w])
for _ in[0]*N:y1,x1,y2,x2=map(int,raw_input().split());print g(s,y1,x1,y2,x2)-g(s1,y1,x1,y2,x2-1)-g(s2,y1,x1,y2-1,x2)

Solution (Python):
f=lambda s,h,w,a,q:s[h+1][w]+s[h][w+1]-s[h][w]+a[h][w]*q
g=lambda s,y1,x1,y2,x2:s[y2][x2]-s[y2][x1-1]-s[y1-1][x2]+s[y1-1][x1-1]
H,W,N=map(int,raw_input().split())
a=[map(int,raw_input()) for _ in [0]*H]
s,s1,s2=[ [[0]*(W+1) for _ in [0]*-~H] for _ in [0]*3 ]
for h in range(H):
	for w in range(W):
		s[h+1][w+1]=f(s,h,w,a,1)
		s1[h+1][w+1]=f(s1,h,w,a,w+1<W and a[h][w+1])
		s2[h+1][w+1]=f(s2,h,w,a,h+1<H and a[h+1][w])
for _ in [0]*N:
	y1,x1,y2,x2=map(int,raw_input().split())
	print g(s,y1,x1,y2,x2)-g(s1,y1,x1,y2,x2-1)-g(s2,y1,x1,y2-1,x2)

Solution (Python):
f=lambda s,y1,x1,y2,x2:s[y2][x2]-s[y2][x1-1]-s[y1-1][x2]+s[y1-1][x1-1]
H,W,N=map(int,raw_input().split())
a=[map(int,raw_input()) for _ in [0]*H]
s,s1,s2=[ [[0]*(W+1) for _ in [0]*-~H] for _ in [0]*3 ]
for h in range(H):
	for w in range(W):
		s[h+1][w+1]=s[h+1][w]+s[h][w+1]-s[h][w]+a[h][w]
		s1[h+1][w+1]=s1[h+1][w]+s1[h][w+1]-s1[h][w]+(a[h][w]*(w+1<W and a[h][w+1]))
		s2[h+1][w+1]=s2[h+1][w]+s2[h][w+1]-s2[h][w]+(a[h][w]*(h+1<H and a[h+1][w]))
for _ in [0]*N:
	y1,x1,y2,x2=map(int,raw_input().split())
	print f(s,y1,x1,y2,x2)-f(s1,y1,x1,y2,x2-1)-f(s2,y1,x1,y2-1,x2)

Solution (Python):
import numpy as np

N, M, Q = map(int, input().split())

S = np.zeros((N+1, M+1), dtype='float32')

for n in range(1, N+1):
    S[n, 1:] = list(map(int, list(input())))


N_cum = np.cumsum(np.cumsum(S, axis=0), axis=1) 

ER_cum = np.cumsum(np.cumsum(S[:N, :]*S[1:, :], axis=0), axis=1) 

EC_cum = np.cumsum(np.cumsum(S[:, :M]*S[:, 1:], axis=0), axis=1) 

sol_list = []

for _ in range(Q):
    x1, y1, x2, y2 = map(int, input().split())
    
    N_nodes = N_cum[x2, y2] - N_cum[x1-1, y2] - N_cum[x2, y1-1] + N_cum[x1-1, y1-1]
    
    N_edges1 = ER_cum[x2-1, y2] - ER_cum[x1-1, y2] - ER_cum[x2-1, y1-1] + ER_cum[x1-1, y1-1]
    
    N_edges2 = EC_cum[x2, y2-1] - EC_cum[x1-1, y2-1] - EC_cum[x2, y1-1] + EC_cum[x1-1, y1-1]
    
    sol_list.append(N_nodes - N_edges1 - N_edges2)

for q in range(Q):
    print(int(sol_list[q]))