text:
## September 17, 2020 Single Round Match 790 Editorials

#### Div2-Easy: Alice’s Birthday

For a given $$K$$, we have to partition the first $$K$$ fibonacci numbers in two sets with equal sum. The fact that $$F_k = F_{k-1} + F_{k-2}$$ means that if $$K$$ is a multiple of 3, then I can give $$F_k$$ to Charlie and $$F_{k-1}, F_{k-2}$$ to Eric, $$F_{k-3}$$ to Charlie and $$F_{k-4}, F_{k-5}$$ to Eric, so on and so forth eventually giving $$F_3$$ to Charlie and $$F_2, F_1$$ to Eric. The other two cases are $$K \equiv 2 \pmod{3}$$ and $$K \equiv 1 \pmod{3}$$.

For the first case, let us try to re-use our previous strategy for the last $$K – 2$$ boxes starting with $$F_k$$ going to Charlie and $$F_{k-1}, F_{k-2}$$ to Eric. In the end, we will be left with the first two boxes containing $$F_1$$ and $$F_2$$, respectively. Since, $$F_1 = F_2 = 1$$ so we can simply give one box to Charlie and other one to Eric.

For the second case, we can try our previous strategy but will always be left with one extra box. The fact that the sum of $$F_1, F_2, \ldots, F_k$$ is odd means that we cannot divide them evenly among Charlie and Eric. Had it been possible to divide it equally, it would mean that the total sum is even, which is not possible and can be proven using induction.

Since we only have to run a loop till $$K$$, so time complexity is $$O(K)$$.

Reference sol in Java

import java.io.*;
import java.util.*;
import javafx.util.Pair;

public class AlicesBirthday {
public int[] partition(int k) {
int[] arr = {-1} ;
if(k % 3 == 1) return arr ;

ArrayList<Integer> charlie = new ArrayList<Integer>(), eric = new ArrayList<Integer>();
while(k > 2) {
}
if(k == 2) {
}
int[] ret = new int[charlie.size()]; int index = 0;
for (Integer value: charlie)
ret[index++] = value;
return ret;
}
}

Reference sol in Python

class AlicesBirthday:
def partition(self, K):
F = [1, 1]
while len(F) < K: F.append( F[-1] + F[-2] )
F = F[:K]
charlie, eric, answer = 0, 0, []
for k in reversed(range(K)):
if charlie > eric:
eric += F[k]
else:
charlie += F[k]
if charlie == eric:
return [-1]

#### Div2-Medium: Bob the Builder

Suppose we are at the given height $$B$$. At any step, we can either move to one of its factors or add $$K$$ to it. While we can move to its factors without any cost, we have to spend $$1$$\$ to add $$K$$. So, basically, we have a given state, few possible options from there, and we want to find out end whether a specific end state (i.e. height $$H$$) is reachable or not. If reachable, we have to do so in minimum cost. This naturally appeals to be modeled by graphs with states representing the nodes and edges representing the transitions. Thus, we have to find the shortest path in a graph, so Dijkstra definitely comes as the first solution in mind.

One may use Dijkstra to solve the problem, but it is not immediately clear when to stop exploring if there does not exist any path. We would never actually need to go beyond $$10^6$$ (proof is shown below). But going so high, Dijkstra may get TLE. The fact that edge weights are only $$0$$ and $$1$$ means that we can use $$0 – 1\ BFS$$ instead of Dijkstra and the time limit was high enough to allow that in case exploration goes so high. However, in practice, the max step explored before reaching a goal seems to be quite low (the worst case we could find was not more than $$5 \cdot 10^4$$); hence normal Dijkstra passes easily.

Reference sol in Java using 0-1 BFS

import java.util.*;
import java.lang.*;
import java.io.*;

public class BobtheBuilder {
public int N = (int) 5e5 + 2 ;
public int maxStep = -1 ;
public List <Integer>[] factors = new ArrayList[N];
public int dist[] = new int[N] ;
public boolean seen[] = new boolean[N] ;

void prep() {
for(int i = 1; i < N; i++) {
factors[i] = new ArrayList<Integer> ();
dist[i] = (int) 2e9 ; seen[i] = false ;
int sq = (int) Math.sqrt(i);
for(int j = 1; j <= sq; j++) {
if(i % j == 0) {
}
}
}
}
public int minimumPrice(int b, int k, int h) {
prep();
Deque<Integer> deq = new ArrayDeque<Integer>(); // seen[b] = true ;
dist[b] = 0 ;
while(deq.size() > 0) {
int top = deq.removeFirst(); if(top > maxStep) maxStep = top ;
if(top == h) return dist[top] ;
for(int f : factors[top]) {
if(dist[top] < dist[f]) {
dist[f] = dist[top] ;
}
}
if( (top + k) < N) {
if( (dist[top] + 1) < dist[top + k]) {
dist[top + k] = dist[top] + 1 ;
}
}
}
return -1 ;
}
}

Reference sol in C++ using Dijkstra

#include <bits/stdc++.h>
using namespace std;

struct record {
int id;
int dist;
};

bool operator<(const record &A, const record &B) {
if (A.dist != B.dist) return A.dist < B.dist;
return A.id < B.id;
}

int solve(int start, int add, int goal) {
if (start % div) return -1;

vector<int> D(1002000,2000000);
D[start] = 0;
set<record> Q;
Q.insert( {start, 0} );
while (!Q.empty()) {
record kde = *Q.begin();
Q.erase( Q.begin() );
if (kde.id >= 1002000) continue;
if (kde.id == goal) break;
for (int div=1; div*div<=kde.id; ++div) if (kde.id % div == 0) {
vector<int> opts = {div, kde.id/div};
for (int kam : opts) {
int ndist = kde.dist;
if (ndist >= D[kam]) continue;
Q.erase( {kam, D[kam]} );
D[kam] = ndist;
Q.insert( {kam, D[kam]} );
}
}
vector<int> opts = { kde.id + add };
for (int kam : opts) {
int ndist = kde.dist + 1;
if (ndist >= D[kam]) continue;
Q.erase( {kam, D[kam]} );
D[kam] = ndist;
Q.insert( {kam, D[kam]} );
}
}
if (D[goal] == 2000000) return -1;
return D[goal];
}

struct BobtheBuilder {
int minimumPrice(int B, int K, int H) { return solve(B,K,H); }
};

Proof
If we are at some state $$C$$, then at anytime, we can go to one of its factors (possibly $$C$$ itself) or add $$K$$. We do this starting from $$B$$ until we reach $$H$$. Suppose it is possible to reach $$H$$, then thinking about the process in reverse is basically from $$H$$ we can go to one of its multiple (possibly $$H$$ itself) or subtract $$K$$ and continue this loop until we reach $$B$$. Thus, we can model this in the following form: $$( .. (((Hx_1 – K)x_2 – K)x_3 – K) ..) = B$$ Simplifying this, we reduce it to a well-known form: $$Hx – Ky = B \quad (\text{for some } x, y \in \mathbb{Z})$$ So, we end up with a linear diophantine equation. Let $$g\ =\ gcd(H, K)$$. This is only solvable when $$g\ |\ B$$ so we can instantly return -1 if $$B$$ is not a multiple of $$g$$. The fact that $$H, K$$ have opposite signs means that this admits infinite positive solutions. From Bezuot’s Identity, we can always find a pair $$s, t$$ such that $$Hs – Kt = g$$ and $$|s| \le K/2\ \&\ |t| \le H/2$$. So we can always find $$s, t$$ lie in the range $$[-500, 500]$$ and since infinite positive solutions are possible and we can shift $$s, t$$, so it’s always possible to find a pair $$(s, t)$$ with both lying in $$[0, 1000]$$. Since $$H, K \le 10^3$$, so the max state we need to explore to ensure we always find a solution is bounded by $$10^6$$. Once we reach $$g$$ from $$H$$, we can repeat the same process to reach $$B$$ from $$g$$. While we did the derivation in the reverse form, it obviously holds in the forward manner starting from $$B$$. Some useful links:

#### Div2-Hard / Div1-Easy: The Social Network

We are given an undirected connected graph $$G$$ with exponential weights of the form $$2^x$$, and we have to find its minimum cut. The problem may look intimidating but what is in favor for us is the fact that all edge weights are distinct.

Let the given weights sorted in decreasing order be $$w_1, w_2, \ldots w_m$$ (i.e. edge weight is $$2^{w_i}$$). Then, for any valid $$j$$ between $$1$$ and $$m – 1$$, we have

So, we would like to keep the edge with the highest weight out of the cut edges set at the cost of ending up with all other edges in the cut since their total sum is strictly less than the highest weight. Thus, we can keep adding edges one by one in the graph (in decreasing order of their weight) as long as one whole connected component does not form.

The edges which do not get added in the graph form our cut edges, and their sum is the answer. This can be done by DSU but was not required since the constraints were significantly small. If we do use DSU, then sorting dominates the runtime with total complexity as $$O(M \text{log}(M))$$.

Reference Sol in Java

import java.util.*;
import java.lang.*;
import java.io.*;

public class TheSocialNetwork {

class Edge implements Comparable<Edge> {
public int u;
public int v;
public int w;

public Edge(int uu, int vv, int cc) {
this.u = uu;
this.v = vv;
this.w = cc;
}

public int compareTo(Edge c) {
if (w < c.w) return 1;
if (w > c.w) return -1;
return 0;
}
}

int M = 1002 ; int C = 100002 ; int N = 302 ;
List <Edge> edges = new ArrayList<Edge>() ;
int[] par = new int[N] ; int[] sz = new int[N] ; int[] pwr = new int[C] ;
int MOD = (int) (1e9 + 7); int comp ;

boolean[] seen = new boolean[C] ;
List < Integer >[] adj = new ArrayList[N] ;
boolean[] visited = new boolean[N] ;
int connected ;

void prepare() {
pwr[0] = 1 ; seen[0] = false ;
for(int i = 1; i < C; i++) {
pwr[i] = (pwr[i-1] * 2) % MOD ;
seen[i] = false ;
}
for(int i = 1; i < N; i++) {
par[i] = i ; sz[i] = 1 ;
visited[i] = false ;
}
}

int getPar(int v) {
if(par[v] == v) return v ;
else return par[v] = getPar(par[v]);
}

int dsu(int u, int v) {
int parU = getPar(u) ; int parV = getPar(v);
if(parU == parV) return 0;

if(sz[u] < sz[v]) { // swap
u ^= v ; v ^= u ; u ^= v ;
}
sz[u] += sz[v] ;
par[parV] = parU ;
return 1 ;
}

public int minimumCut(int n, int m, int[] u, int[] v, int[] l) {
prepare();
for(int i = 0; i < m; i++) {
Edge e = new Edge(u[i], v[i], l[i]);
}
Collections.sort(edges);
int ans = 0; comp = n ;
for(Edge e : edges) {
System.out.println(e.w);
int parU = getPar(e.u) ; int parV = getPar(e.v);
if(parU == parV) {
}
else {
if(comp > 2) {
comp-- ;
dsu(e.u, e.v);
}
else {
ans += pwr[e.w] ;
ans %= MOD ;
}
}
}
return ans ;
}
}

Reference Sol in Python

def pow2(n):
if n == 0: return 1
t = pow2(n//2)
t *= t
if n % 2: t *= 2
return t % 1000000007

class TheSocialNetwork:
def minimumCut(self, n, m, u, v, l):
edges = [ (-z,x-1,y-1) for x,y,z in zip(u,v,l) ]
edges.sort()
component_count = n
component = list(range(n))
for w,x,y in edges:
if component[x] == component[y]: continue
if component_count == 2:
else:
cx, cy = component[x], component[y]
component = [ cx if c==cy else c for c in component ]
component_count -= 1

Side Thought: We were wondering whether it is possible to solve this problem using max-flow?

#### Div1-Medium: Proposal Optimization

Such kind of ratio problems can be solved by using a niche trick. I think it is hard to develop this line of thinking unless you have encountered it before.

Let us say the optimal ratio is $$Q$$. Thus, $$\frac{\sum R_i}{\sum T_i} = Q$$. What we can do is a binary search on $$Q$$, and check if a ratio is feasible or not. So, for a given $$Q$$, we have to check whether $$\frac{\sum R_i}{\sum T_i} \ge Q$$ is possible, which is the same thing as $$\sum R_i – Q*\sum T_i \ge 0$$ (let’s call the L.H.S. $$P$$). The constraints on $$R_i, T_i$$ seem to be quite high to do any form of *DP* while constraints on the grid’s dimension, $$N, M$$ are quite low although obviously not permitting brute force solutions. However, meet-in-the-middle looks promising since we can store $$P$$ while doing brute force from one corner of the grid and checking from another (meeting in cells lying on some sort of diagonal of the grid).

In fact, it is the intended solution, and we also have to take care of the costs along with $$P$$ to actually know if any of the options is actually feasible. So, when doing brute from one side, we store $$(P_i, Cost_i)$$ pairs on those diagonal cells and later sort them according to $$Cost_i$$ to ease our process in the next brute. While doing brute from the other end, we have $$(P_j, Cost_j)$$, and we want the highest $$P_i$$ whose cost satisfies $$Cost_i \le K – Cost_j$$. By storing $$P_i$$ in prefix maximum manner, we can binary search on the highest index option with cost $$\le K – Cost_j$$ and take its $$P_i$$ (since we already sorted) to get our answer.

The best diagonal cells for a grid of size $$(N, M)$$ would be the set $${(i, j) : i + j == (N+M)/2 }$$. Since $$N$$x$$M \le 300$$, so the number of such diagonal cells (let us denote that by $$C$$) can’t exceed $$17 \quad (\sqrt{300} = 17 )$$. We can focus on the analysis of square matrix since it will dominate the rectangular grids. Let the number of binary search iterations required to find the optimal ratio be $$L$$ (80 iterations suffice for $$10^{-9}$$ error). To reach a diagonal cell $$(i, j)$$ there are $$i + j \choose i$$ ways and then we have to binary search on those many paths which we got from first brute force. Thus, doing meet-in-the-middle takes $$\sum^C_{i = 1} {C \choose i} \cdot log( {C \choose i} )\ \le \ C \cdot 2^C$$ Hence, the total time required is bounded by $$O(L \cdot C \cdot 2^C)$$.

Reference Sol in C++

#include <bits/stdc++.h>
#include<ext/pb_ds/assoc_container.hpp>
#include<ext/pb_ds/tree_policy.hpp>

using namespace std ;
using namespace __gnu_pbds;

template <typename T> // *s.find_by_order(0), s.order_of_key(2) ;
using ordered_set = tree<T, null_type, less<T>, rb_tree_tag, tree_order_statistics_node_update>;

#define reMin(a, b) a = min(a, b)
#define reMax(a, b) a = max(a, b)

#define lint long long
#define pb push_back
#define F first
#define S second
#define sz(x) (int)x.size()
#define all(x) begin(x), end(x)
#define SET(x, val) memset(x, val, sizeof(x))
#define fastio ios_base::sync_with_stdio(0); cin.tie(0); cout.tie(0)

typedef vector < int > vi ;
typedef pair < int, int > pii ;

const int N = 300 + 2 ;
const int MOD = 1e9 + 7 ;
const lint INF = 1e18 ;

double best = 0 ;
int n, m, k ;
int roses[N][N], tulips[N][N], cost[N][N] ;

vector < pair < int, double > > meet[N][N] ;
vector < double > ratios[N][N] ;
int diagonal ;
bool possible = 0 ;
bool ratioSatisfied = 0 ;

void down(int i, int j, int rose, int tulip, int totalCost, double& R) {
if(i >= n or j >= m) return ;
rose += roses[i][j] ; tulip += tulips[i][j] ; totalCost += cost[i][j] ;
if(totalCost > k) return ;
if(i + j == diagonal) {
double cur = rose ; cur -= R * tulip ;
meet[i][j].pb({totalCost, cur});
return ;
}
down(i+1, j, rose, tulip, totalCost, R); down(i, j+1, rose, tulip, totalCost, R);
}

void up(int i, int j, int rose, int tulip, int totalCost, double& R) {
if(i < 0 or j < 0) return ;
rose += roses[i][j] ; tulip += tulips[i][j] ; totalCost += cost[i][j] ;
if(totalCost > k) return ;
if(i + j == diagonal) {
rose -= roses[i][j] ; tulip -= tulips[i][j] ; totalCost -= cost[i][j] ;
double cur = rose ; cur -= R * tulip ;

pair < int, double > item = {k - totalCost, 2e9} ;
int idx = upper_bound(all(meet[i][j]), item) - meet[i][j].begin();
idx-- ;
if(idx >= 0) {
possible = 1 ; assert(idx < sz(ratios[i][j]) and sz(ratios[i][j]) == sz(meet[i][j]));
double mx = ratios[i][j][idx] ;

if((mx + cur) >= 0) // should be EPS!!
ratioSatisfied = 1 ;
}
return ;
}
up(i-1, j, rose, tulip, totalCost, R); up(i, j-1, rose, tulip, totalCost, R);
}

void check(double ratio) {
down(0, 0, 0, 0, 0, ratio);

for(int i = 0; i < n; i++) {
int j = diagonal - i ;
if(j >= m or j < 0) continue ;

if(sz(meet[i][j]) <= 0) continue ;
sort(all(meet[i][j]));

ratios[i][j].pb(meet[i][j][0].S);
for(int p = 1; p < sz(meet[i][j]); p++)
ratios[i][j].pb(max(meet[i][j][p].S, ratios[i][j].back()));
}
up(n-1, m-1, 0, 0, 0, ratio);
}

double solve() {
double low = 0, high = 1e6 ; ratioSatisfied = 0 ;

for(int iter = 0; iter < 70; iter++) {
double mid = (low + high) / 2 ;
check(mid);
if(ratioSatisfied) low = mid ;
else high = mid ;

// clear stuff
ratioSatisfied = 0 ;
for(int i = 0; i < n; i++) {
int j = diagonal - i ;
if(j < 0 or j >= m) continue ;
ratios[i][j].clear(); meet[i][j].clear();
}
}
return low ;
}

class ProposalOptimization {
public:
double bestPath(int rows, int cols, int K, vi roz,  vi tolip, vi expense);
};

double ProposalOptimization::bestPath(int rows, int cols, int K, vi roz,  vi tolip, vi expense ){
n = rows ; m = cols ; k = K ;
for(int i = 0; i < n; i++) for(int j = 0; j < m; j++) roses[i][j] = roz[i*cols + j] ;
for(int i = 0; i < n; i++) for(int j = 0; j < m; j++) tulips[i][j] = tolip[i* cols + j] ;
for(int i = 0; i < n; i++) for(int j = 0; j < m; j++) cost[i][j] = expense[i*cols + j] ;
double best = solve();
if(!possible) return -1 ;
else return best ;
}

Reference Sol in Java

import java.util.*;
import java.math.*;

public class ProposalOptimization {
class Option implements Comparable<Option> {
public int orchids;
public int tulips;
public int cost;

public Option(int orchids, int tulips, int cost) {
this.orchids = orchids;
this.tulips = tulips;
this.cost = cost;
}

public int compareTo(Option c) {
if (cost < c.cost) return -1;
if (cost > c.cost) return 1;
if (orchids < c.orchids) return -1;
if (orchids > c.orchids) return 1;
if (tulips < c.tulips) return -1;
if (tulips > c.tulips) return 1;
return 0;
}
}

Option[][] getOptions(int[][] orchids, int[][] tulips, int[][] costs, int steps, boolean secondPass) {
int R = orchids.length, C = orchids[0].length;
int[][] options = new int[][] { {0,0,0,0,0} };
for (int d=1; d<=steps; ++d) {
int newOptionCount = 0;
for (int i=0; i<options.length; ++i) {
if (options[i][0] + 1 < R) ++newOptionCount;
if (options[i][1] + 1 < C) ++newOptionCount;
}
int[][] newOptions = new int[newOptionCount][5];
for (int i=0, j=0; i<options.length; ++i) {
if (options[i][0] + 1 < R) {
int r = newOptions[j][0] = options[i][0] + 1;
int c = newOptions[j][1] = options[i][1];
int mul = (secondPass && d == steps) ? 0 : 1;
newOptions[j][2] = options[i][2] + mul * orchids[r][c];
newOptions[j][3] = options[i][3] + mul * tulips[r][c];
newOptions[j][4] = options[i][4] + mul * costs[r][c];
++j;
}
if (options[i][1] + 1 < C) {
int r = newOptions[j][0] = options[i][0];
int c = newOptions[j][1] = options[i][1] + 1;
int mul = (secondPass && d == steps) ? 0 : 1;
newOptions[j][2] = options[i][2] + mul * orchids[r][c];
newOptions[j][3] = options[i][3] + mul * tulips[r][c];
newOptions[j][4] = options[i][4] + mul * costs[r][c];
++j;
}
}
options = newOptions;
}
for (int i=0; i<options.length; ++i) ++answerCounts[options[i][0]];
for (int r=0; r<R; ++r) answer[secondPass ? R-1-r : r] = new Option[answerCounts[r]];
for (int i=0; i<options.length; ++i) {
int r = options[i][0];
answer[secondPass ? R-1-r : r][ --answerCounts[r] ] = new Option( options[i][2], options[i][3], options[i][4] );
// System.out.println("row "+r+" added option "+options[i][2]+" "+options[i][3]+" "+options[i][4] );
}
}

int[][] flip(int[][] array) {
int R = array.length, C = array[0].length;
for (int r=0; r<R; ++r) for (int c=0; c<C; ++c) answer[r][c] = array[R-1-r][C-1-c];
}

boolean solvable(double x, Option[] options1, Option[] options2, int K) {
int b = 0;
double bestb = -1e20;
for (int a=options1.length-1; a>=0; --a) {
while (b < options2.length && options1[a].cost + options2[b].cost <= K) {
bestb = Math.max( bestb, options2[b].orchids - x * options2[b].tulips );
++b;
}
if (options1[a].orchids - x * options1[a].tulips + bestb >= 0) return true;
}
return false;
}

public double bestPath(int R, int C, int K, int[] _orchids, int[] _tulips, int[] _costs) {

int[][] orchids = new int[R][C];
int[][] tulips = new int[R][C];
int[][] costs = new int[R][C];
for (int r=0; r<R; ++r) for (int c=0; c<C; ++c) {
orchids[r][c] = _orchids[r*C + c];
tulips[r][c] = _tulips[r*C + c];
costs[r][c] = _costs[r*C + c];
}

Option[][] options1 = getOptions( orchids, tulips, costs, (R+C-2)/2, false);
Option[][] options2 = getOptions( flip(orchids), flip(tulips), flip(costs), (R+C-1)/2, true);

for (int i=0; i<options1.length; ++i) Arrays.sort(options1[i]);
for (int i=0; i<options2.length; ++i) Arrays.sort(options2[i]);

int bestCost = K+1;
for (int i=0; i<options1.length; ++i) if (options1[i].length > 0 && options2[i].length > 0) bestCost = Math.min( bestCost, options1[i][0].cost + options2[i][0].cost );
if (bestCost > K) return -1;

double lo = 0, hi = 1e6 + 1;
for (int step=0; step<200; ++step) {
double med = (lo+hi) / 2;
boolean isSolvable = false;
for (int i=0; i<options1.length; ++i) if (solvable(med,options1[i],options2[i],K)) {
isSolvable = true;
break;
}
if (isSolvable) lo = med; else hi = med;
}
return lo;
}
}

#### Div1-Hard: Tale Of Two Squares:

A number can be expressed as a sum of two sqaures if and only if all its prime factors of the form $$4k+3$$ occur even times. So while taking the product of few numbers, we are only interested in its prime decomposition; hence, we can store a bit-vector for each $$A_i$$ denoting odd/even occurrence of such prime factors.

From now on, whenever I talk about primes it is assumed to be of the form $$4k + 3$$ (for some $$k \in \mathbb{Z}$$) and that the reader is familar with basis of a vector space and gaussian elimination. Please go through this if you are not familiar with the usage of those concepts in CP:
Introductory Blog on using Linear Algebraic techniques to solve XOR related problems

Now, $$A_i \le 10^7$$ so #primes is around $$3* 10^5$$. Normal Gaussian elimination would easily TLE, even with bitsets. Still, lets imagine doing Gaussian elimination traversing from the highest bit to the lowest and xoring when necessary. Notice that, when you are traversing, if any prime factor > $$\sqrt A_i$$ exists (at max one can exist), the moment you xor it, no more bits > $$\sqrt A_i$$ will be on anymore (since only one such was on earlier and now you have xored it). Now, you can simply simulate Gaussian elimination on lower-order primes (i.e. $$primes \le \sqrt A_i$$) and handle that one highest bit (if it exists) separately. This way, combined with bitset for lower-order primes, solves the problem without queries.

For queries, let us process them offline. Let $$queries[l]$$ store all the queries with their left endpoint as $$l$$. We will insert $$A_i$$ into our basis from the right side, $$i$$ going from $$N-1$$ to $$0$$. After we have inserted the element at index $$i$$, we will answer all the queries which begin their, i.e., the ones in $$queries[i]$$.

We will use a simple fact about linear combination of vectors but in a clever manner to answer the queries. Suppose we are about to try and insert $$v$$ into our $$Basis$$ (may not insert if it turns out to be dependant). If $$v = b_1 + b_2 + \cdots b_k$$ (where $$b_i \in Basis$$), then we can remove one of the $$b_i$$ and instead insert $$v$$ into its place, keeping the basis intact. So, we will use this fact while inserting — we will try to keep such elements in the basis which have leftmost index so as to answer as many queries as possible. We can do this greedily: while inserting $$v$$, whenever we are about to xor it with some $$b_i$$, we will first check which has the least index and $$swap(b_i, v)$$ if $$v$$‘s index is smaller than $$b_i$$. If we do end up swapping, then $$v$$ becomes a part of the basis (in place of $$b_i$$), and we progress doing Gaussian Elimination with $$b_i$$ instead. Otherwise, we move forward with $$v$$ like we normally do. Please refer to the attached code for more details.

Now, when we answer queries which have left index $$i$$, we will have the basis vectors from $$A[i:N]$$ having the least possible index. So for a query $$i\ R$$, we simply need to count how many vectors in our basis have index $$\le R$$. This can be easily done with some BIT or set.

Let $$K$$ denote the number of primes $$\le \sqrt A_i \ ( = 225)$$. We store a bitset of size $$K$$ for every prime. While doing Gaussian elimination, we handle the large bit manually while following the normal procedure for lower-order primes; hence it takes $$O(K^2/64)$$ time to insert a $$A_i$$. After that, we answer queries at index $$i$$ in $$O(log(N))$$ time by querying in our BIT. Hence, the overall complexity is $$O(N\cdot \frac{K^2}{64} + Q\cdot log(N))$$, dominated by Gaussian elimination using bitsets:

Reference Sol in Java

import java.util.*;
import java.lang.*;
import java.io.*;

public class TaleOfTwoSquares {

class Query {
public int L;
public int R;
public int idx;

public Query(int l, int r, int j) {
this.L = l;
this.R = r;
this.idx = j;
}
}

int K = (int) 1e7 + 2 ;  int N = 100002 ; int Q = 100002 ;
int BS = 256 ; int SQRT = 3200 ; int MAX_PRIMES = (int) 4e5 ;

int spf[] = new int[K] ;
List<Integer> arr = new ArrayList<Integer>();

Map<Integer, Integer> map = new HashMap<Integer, Integer>();
Map<Integer, Integer> highMap = new HashMap<Integer, Integer>();

boolean filled[] = new boolean[SQRT] ;

int indices[] = new int[K] ;
List<BitSet> basis = new ArrayList<BitSet>();
List<BitSet> basisHigh = new ArrayList<BitSet>();

List<Query>[] queries = new ArrayList[Q] ;
int basisIndex[] = new int[K] ;
int bit[] = new int[N] ;
int pwr2[] = new int[MAX_PRIMES] ;

int MOD = 998244353 ;
Scanner s = new Scanner(System.in);

void spfsieve()
{
for(int i = 2; i < K; i++) spf[i]=-1;
int cnt = 0 ;
for(int i = 2; i < K; i++)
{
if(spf[i] != -1) continue;
if(i % 4 == 3) indices[i] = cnt++ ;
for(int j = i; j < K; j+=i)
if(spf[j]==-1)
spf[j]=i;
}
assert cnt < MAX_PRIMES ;
}

void update(int i, int x) {
for(; i < N; i += i & -i)
bit[i] += x ;
}

int getPrefixSum(int i) {
int s = 0;
for(; i > 0; i -= i & -i)
s += bit[i] ;
return s ;
}

void prepare() {
for(int i = 0; i < Q; i++)
queries[i] = new ArrayList < Query > ();
pwr2[0] = 1 ;
for(int i = 0; i < MAX_PRIMES; i++) {
BitSet lst = new BitSet(BS);
if(i > 0)
pwr2[i] = (pwr2[i-1] * 2) % MOD ;
}
for(int i = 0; i < SQRT; i++) {
filled[i] = false ;
BitSet lst = new BitSet(BS);
}
}

int insertVector(int highest, int idx, BitSet bs) {
// handle highest bit manually
if(highest >= SQRT) {
int bitIndex = indices[highest] ;

if(highMap.containsKey(bitIndex)) {

if(idx < basisIndex[bitIndex]) {
update(basisIndex[bitIndex], -1); update(idx, 1);

// swap
bs.xor(basisHigh.get(bitIndex));
BitSet cpy = basisHigh.get(bitIndex);
cpy.xor(bs); basisHigh.set(bitIndex, cpy);
bs.xor(basisHigh.get(bitIndex));

int tmp = basisIndex[bitIndex] ;
basisIndex[bitIndex] = idx ;
idx = tmp;
}
bs.xor(basisHigh.get(bitIndex));
}
else {
highMap.put(bitIndex, idx);
basisHigh.set(bitIndex, bs);
basisIndex[bitIndex] = idx ; update(idx, 1);
return 0;
}
}
// gauss
for(int i = SQRT - 1; i >= 0; i--) {
if(!bs.get(i)) continue ;
if(!filled[i]) {
filled[i] = true ;
basis.set(i, bs) ; update(idx, 1);
basisIndex[i] = idx ;
return 0 ;
}
if(idx < basisIndex[i]) {
update(basisIndex[i], -1); update(idx, 1);

// swap
bs.xor(basis.get(i));
BitSet cpy = basis.get(i); cpy.xor(bs); basis.set(i, cpy);
bs.xor(basis.get(i));

int tmp = basisIndex[i] ;
basisIndex[i] = idx ;
idx = tmp;
}
bs.xor(basis.get(i));
}
return 1;
}

public int count(int n, int[] aprefix, int query_cnt, int[] ql, int[] qr, int seed) {
spfsieve(); prepare(); long state = seed ; long modo = 1 ; modo <<= 31 ;

for(int i = aprefix.length; i < n; i++) {
state = (state * 1103515245 + 12345) % modo ;
arr.add( (int) (1 + (state % 10000000)) );
}

int qprefix = ql.length; assert qprefix == qr.length ;
for(int i = 0; i < qprefix; i++) {
Query q = new Query(ql[i] + 1, qr[i] + 1, i);
}

for(int i = qprefix; i < query_cnt; i++) {
state = (state * 1103515245 + 12345) % modo ;
int x = (int) state % n ; x++ ;
state = (state * 1103515245 + 12345) % modo ;
int y = (int) state % n ; y++ ;
int l = Math.min(x, y); int r = Math.max(x, y);
Query q = new Query(l, r, i);
}

int sum = 0 ;
for(int i = n; i > 0; i--) {
int now = arr.get(i - 1) ;
BitSet bs = new BitSet(BS);
int highest = -1 ;

while(now > 1) {
int primeFactor = spf[now] ;
int occurs = 0 ;
while(now % primeFactor == 0) {
now /= primeFactor ;
occurs++ ;
}
if(primeFactor % 4 == 3 && occurs % 2 == 1) {
// System.out.println(arr.get(i));
if(primeFactor >= SQRT) {
assert highest == -1 ;
highest = primeFactor ;
}
else {
int pos = indices[primeFactor] ;
bs.set(pos);
}
}

}
insertVector(highest, i, bs);
for(Query qq : queries[i]) {
int sz = getPrefixSum(qq.R);
int tot = qq.R - qq.L + 1 ;
int dependant = tot - sz ;
// System.out.println(qq.L + " " + qq.R + " " + sz + " - " + dependant);
int ans = pwr2[dependant] - 1 ;
if(ans < 0) ans += MOD ;
sum += ans ; sum %= MOD ;
}
}
return sum ;
}
}

Reference sol in C++

#include <bits/stdc++.h>
using namespace std;

const int SQRTMAX = 3162;
const int SMALL_PRIME_COUNT = 225;
const int SMALL_PRIMES[] = {3, 7, 11, 19, 23, 31, 43, 47, 59, 67, 71, 79, 83, 103, 107, 127, 131, 139, 151, 163, 167, 179, 191, 199, 211, 223, 227, 239, 251, 263, 271, 283, 307, 311, 331, 347, 359, 367, 379, 383, 419, 431, 439, 443, 463, 467, 479, 487, 491, 499, 503, 523, 547, 563, 571, 587, 599, 607, 619, 631, 643, 647, 659, 683, 691, 719, 727, 739, 743, 751, 787, 811, 823, 827, 839, 859, 863, 883, 887, 907, 911, 919, 947, 967, 971, 983, 991, 1019, 1031, 1039, 1051, 1063, 1087, 1091, 1103, 1123, 1151, 1163, 1171, 1187, 1223, 1231, 1259, 1279, 1283, 1291, 1303, 1307, 1319, 1327, 1367, 1399, 1423, 1427, 1439, 1447, 1451, 1459, 1471, 1483, 1487, 1499, 1511, 1523, 1531, 1543, 1559, 1567, 1571, 1579, 1583, 1607, 1619, 1627, 1663, 1667, 1699, 1723, 1747, 1759, 1783, 1787, 1811, 1823, 1831, 1847, 1867, 1871, 1879, 1907, 1931, 1951, 1979, 1987, 1999, 2003, 2011, 2027, 2039, 2063, 2083, 2087, 2099, 2111, 2131, 2143, 2179, 2203, 2207, 2239, 2243, 2251, 2267, 2287, 2311, 2339, 2347, 2351, 2371, 2383, 2399, 2411, 2423, 2447, 2459, 2467, 2503, 2531, 2539, 2543, 2551, 2579, 2591, 2647, 2659, 2663, 2671, 2683, 2687, 2699, 2707, 2711, 2719, 2731, 2767, 2791, 2803, 2819, 2843, 2851, 2879, 2887, 2903, 2927, 2939, 2963, 2971, 2999, 3011, 3019, 3023, 3067, 3079, 3083, 3119 };
vector<int> SMALL_PRIME_ID, pow2;

struct Fenwick1D { // {{{
int size;
vector<int> T;

Fenwick1D(int maxval) {
size = 1;
while (size < maxval) size <<= 1;
T.clear();
T.resize(size+1,0);
}

void update(int x, int delta) { // assumes 1 <= x <= init_maxval
while (x <= size) { T[x] += delta; x += x & -x; }
}

int sum(int x1, int x2) { // sum in the closed interval [x1,x2]
int res=0;
--x1;
while (x2) { res += T[x2]; x2 -= x2 & -x2; }
while (x1) { res -= T[x1]; x1 -= x1 & -x1; }
return res;
}

int find(int sum) { // largest z such that sum( [1,z] ) <= sum
int idx = 0, bitMask = size;
while (bitMask && (idx < size)) {
int tIdx = idx + bitMask;
if (sum >= T[tIdx]) { idx=tIdx; sum -= T[tIdx]; }
}
return idx;
}
}; // }}}

void init() {
SMALL_PRIME_ID.clear();
SMALL_PRIME_ID.resize(SQRTMAX+1, -1);
for (int i=0; i<SMALL_PRIME_COUNT; ++i) SMALL_PRIME_ID[ SMALL_PRIMES[i] ] = i;
pow2.clear();
pow2.push_back(1);
for (int i=1; i<100005; ++i) pow2.push_back( (pow2.back()*2) % 998244353 );
}

struct signature {
bitset< SMALL_PRIME_COUNT > small_bits;
int large_bit;
signature() : large_bit(-1) {}

int highest_prime() {
if (large_bit != -1) return large_bit;
for (int i=SMALL_PRIME_COUNT-1; i>=0; --i) if (small_bits.test(i)) return SMALL_PRIMES[i];
return 0;
}

void do_xor(const signature &other) {
if (other.large_bit != -1) {
if (large_bit != -1) {
assert( large_bit == other.large_bit );
large_bit = -1;
} else {
large_bit = other.large_bit;
}
}
small_bits ^= other.small_bits;
}
};

ostream& operator<< (ostream& out, const signature &S) {
/*
int max = SMALL_PRIME_COUNT - 1;
while (max > 0 && !S.small_bits.test(max)) --max;
for (int i=0; i<=max; ++i) out << int( S.small_bits.test(i) );
*/
if (S.large_bit != -1) out << S.large_bit << "+";
for (int i=9; i>=0; --i) out << int(S.small_bits.test(i));
return out;
}

struct query {
int l, r;
};

bool operator< (const query &A, const query &B) {
if (A.l != B.l) return A.l > B.l;
return A.r < B.r;
}

signature get_signature(int N) {
for (int d=2; d*d<=N; ++d) {
if (N % d) continue;
int cnt = 0;
while (N % d == 0) { ++cnt; N /= d; }
if (d % 4 != 3) continue;
if (cnt % 2 == 0) continue;
}
if (N > 1 && N % 4 == 3) {
}
}

struct TaleOfTwoSquares {
int count(int N, vector<int> Aprefix, int Q, vector<int> Lprefix, vector<int> Rprefix, int seed) {
init();

long long state = seed;

vector<int> A = Aprefix;
while (int(A.size()) < N) {
state = (state * 1103515245 + 12345) % (1LL << 31);
A.push_back(1 + (state % 10000000));
}

if (N <= 20) for (int n=0; n<N; ++n) cout << A[n] << ", "; cout << endl;

vector<int> L = Lprefix, R = Rprefix;
while (int(L.size()) < Q) {
state = (state * 1103515245 + 12345) % (1LL << 31);
int x = state % N;
state = (state * 1103515245 + 12345) % (1LL << 31);
int y = state % N;
L.push_back( min(x,y) );
R.push_back( max(x,y) );
}

if (Q <= 20) for (int n=0; n<Q; ++n) cout << L[n] << ", "; cout << endl;
if (Q <= 20) for (int n=0; n<Q; ++n) cout << R[n] << ", "; cout << endl;

vector<signature> B;
for (int n=0; n<N; ++n) B.push_back( get_signature(A[n]) );

vector< vector<int> > query_by_L(N);
for (int q=0; q<Q; ++q) query_by_L[ L[q] ].push_back( R[q] );

Fenwick1D InBase(N+2);
unordered_map<int,int> base_by_highest_prime;

for (int n=N-1; n>=0; --n) {
InBase.update(N+1,+1);
InBase.update(N+1,-1);
// add number A[n] to the base
if (B[n].highest_prime() != 0) {
InBase.update(n+1,+1);
int where = n;
while (true) {
int p = B[where].highest_prime();
if (p == 0) {
InBase.update(where+1,-1);
break;
}
if (!base_by_highest_prime.count(p)) {
base_by_highest_prime[p] = where;
break;
}
int nxt = base_by_highest_prime[p];
if (nxt > where) {
base_by_highest_prime[p] = where;
B[nxt].do_xor( B[where] );
where = nxt;
} else {
B[where].do_xor( B[nxt] );
}
}
}
// answer queries that begin here
for (int r : query_by_L[n]) {
int length = r-n+1;
int base_size = InBase.sum(n+1,r+1);
}
}

}
};

nikhil_chandak

Guest Blogger

categories & Tags

UNLEASH THE GIG ECONOMY. START A PROJECT OR TALK TO SALES
Close


text:
# Help understanding the proof that a Riemann Integrable function is bounded

Suppose $f$ is Riemann integrable over $[a,b]$. The proof I am reading that shows that $f$ must be bounded, just says if $f$ is unbounded then there is some point $x$, such that for any $n\in\mathbb{N}$ and $\epsilon>0$ if we take any partition divded into $n$ parts we have $|f(x)|>n\epsilon$ and so if we take $x$ as a tag in the partition then $$|S(f,\mathcal{P})-L|>\epsilon\space\space\space\space\space\space (\dagger),$$ where $S(f,\mathcal{P})$ is the Riemann sum of $f$ with respect to the partition $\mathcal{P}$ and $L=\int_a^bf$.

My trouble is trying to get $(\dagger)$, I am not sure how split up the Riemann sum to get the inequality, I was thinking the reverse triangle but get nowhere. So any help will be really appreciated and needed.

• This is direct from the definition (one could even say that boundedness is a part of the definition of Riemann integrability). To wit, if $f$ is Riemann integrable, then $f$ is bounded from above and from below by some step functions. Every step function is bounded hence $f$ is bounded. The "proof" in your post seems unnecessarily convoluted (in addition, it is frankly suspect). – Did Apr 5 '17 at 12:23

You have not told us what your working definition of Riemann integrability is.

If integrability is defined in terms of lower and upper sums, and $f$ is unbounded above, then every single upper sum is undefined (or $=\infty$, if you prefer). Basta; no further epsilontics needed.

If integrability is defined in terms of quantities like $\|\Delta f\|_{I_j}:=\sup_{x, \ y\in I_j}|f(x)-f(y)|$ then again for every partition at least one of these quantities is $=\infty$, hence you can never make $\sum_j\|\Delta f\|_{I_j}\>|I_j|<\epsilon$. But this is required before you can even think of computing "general Riemann sums" $\sum_{j=1}^N f(\tau_j)\>|I_j|$ as approximations to the intended integral.

• My suspicion is that a 'tagged partition' is a sequence $a = t_0< t_1 < ... < t_n = b$ together with values $x_i \in (t_{i-1}, t_i)$, and there's some notion of the sum associated to a tagged partition, a partial order on tagged partitions, and a limit as the mesh of the partition goes to zero. (If this is right, I think it's probably pedagogically a bad way to go...even proving that $f(x) = x$ is integrable on $[0, 1]$ could take a week of class time!) – John Hughes Apr 5 '17 at 14:36
• Yes, that is right, we do however prove that this "classical" definition is equivalent to the Darboux definition where $f$ is taken to be bounded, I do have a proof (not proof by contradiction) that $f$ must be bounded in E.D. Bloch's Real numbers and real analysis bit it is a bit long winded so I found this proof which seems much shorter but my trouble is showing $(\dagger)$. So this proof is from the Real Analysis wikibook (en.wikibooks.org/wiki/Real_Analysis/Riemann_integration) so I apologise for not being clear in my 1st post. – user152874 Apr 5 '17 at 15:10
• So any help with showing $(\dagger)$ will really be appreciated. I am guessing it is probably quiet simple and easier statement are given more detail in the text. But I just can't see it. Thanks again – user152874 Apr 5 '17 at 15:12

This seems wrong to me. I might be misunderstanding the use of "tag", however.

Look at $$f(x) = \frac{1}{x-1} + \frac{1}{x+1}$$ on the interval $(-1, 1)$.

The points $+1$ and $-1$ are examples of the point $x$ with the necessary property (unboundedness near $x$). But if you take the $k$-interval equipartition with its $k+1$ breakpoints evenly spaced from $-1$ to $1$, and pick $k$ equispaced points between $-1$ and $1$ as the "samples" or "tags", you find that $S(f, P) = 0$. That's at least consistent with $L = 0$ (even through the function is not, in fact, Riemann integrable.

In short: the claim that unboundedness messes things up is correct. The alleged proof of this claim appears wrong.

• This the proof from the text: If $f$ is unbounded. For every $n$ divide the interval into $n$ parts. Hence, for every $n$, $f$ is unbounded on at least one of these $n$ parts. Call it $I_n$ Now, let $\epsilon >0$ be given. Consider an arbitrary $\delta >0$. Let $\mathcal{P}$ be a tagged partition such that $\|\mathcal {P}\|<\delta$ and $(I_n,t_n)\in \mathcal{P}$, where $t_n$ is taken so as to satisfy $|f(t_n)|>n\epsilon$. Thus we have that $|S(f,\mathcal{P})-L|>\epsilon$. A contradiction to the fact that $f$ is Riemann integrable. So what are the details to show $(\dagger)$? Thanks – user152874 Apr 5 '17 at 13:11
• What does $\|P\|$ denote? The "mesh" of $P$, i.e., the length of the largest interval in $P$? And what is $L$? It seems to be unbound in this quotation at least. Is it the purported value of the integral? The sentence "Thus we have that ..." seems to hide a great deal (and continuing to look at my example, I don't actually believe it's true). – John Hughes Apr 5 '17 at 13:51
• For $f$ unbounded, what's the definition of the Upper Sum of $f$ with respect to a partition $P$? In my favorite definition, we let $M_i$ be the max of $f$ on the $i$th subinterval of $P$, and the upper sum involves $M_i (t_{i+1} - t_i)$, but that only makes sense because the definition restricts to bounded functions at the start. – John Hughes Apr 5 '17 at 14:01
• A tagged partition of $[a,b]$ is defined as $\{([x_{k-1},x_k]),t_k)\}_0^n$, where $a=x_0<...<x_n=b$ and the "tags" $t_k\in[x_{k-1},x_k]$, with $\mathbb{P}_{[a,b]}$ is the set of all tagged partitions over $[a,b]$. Yes $\|\mathcal{P}\|$ is the mesh of the partition. The Riemann sum of $f$ over $[a,b]$ w.r.t $\mathcal{P}$ is defined as $S(f,\mathcal{P}):=\sum_1^n f(t_k)(x_k-x_{k-1})$ an $f$ is Riemann integrable with $\int_a^bf=L$ iff $$(\forall\epsilon>0)(\exists\delta>0)(\forall \mathcal{P}\in\mathbb{P}_{[a,b]})(\|\mathcal{P}\|<\delta\Rightarrow |S(f,\mathcal{P})-L|<\epsilon)$$ – user152874 Apr 5 '17 at 15:23
• In that case, the proof is certainly wrong, for $I_n$ has length $(b-a)/n$, and for $\delta = (b-a)/(2n)$, it's impossible for $\|P\| = \delta$ and for $I_n \in P$. I'm growing to have grave doubts about your text... – John Hughes Apr 5 '17 at 16:30

In comments, you wrote the following (which I've slightly polished):

A tagged partition of $[a,b]$ is a sequence
$\{([x_{k−1},x_{k}]),t_k)\}_{k=0}^n$, where $a=x_0<...<x_n=b$ and the "tags" $t_k$ satisfy $$t_k\in [x_{k−1},x_k]$$ for $k = 1, \ldots n$.

$\Bbb P[a,b]$ denotes the set of all tagged partitions over $[a,b]$.

$\|P\| = \max_{k=1}^n (x_k - x_{k-1})$ is called the mesh of the partition.

The Riemann sum of $f$ over $[a,b]$ w.r.t a partition $P$ is defined as $$S(f,P):=\sum_1^n f(t_k)(x_k−x_{k−1}).$$

Finally, $f$ is Riemann integrable with $\int_a^b f = L$ iff $$\forall \epsilon > 0~ \exists \delta > 0 ~\forall P \in \Bbb P[a,b], \|P\| < \delta \implies |S(f,P)−L|<\epsilon.$$

You want to prove that if $f$ is unbounded, $f$ is not integrable.

Proof: By contradiction: Suppose $f$ is unbounded, but that it's also integrable with integral $L$.

Pick $\epsilon = 1$. Then there's a $\delta$ with the property that for any tagged partition $P$ with mesh less than $\delta$, $|S(f, P) - L | < 1$. Let $P$, with split points $x_i$ and tags $t_i$, be such a partition** . Let $n$ be the number of tags in $P$.

I'll exhibit a different partition, $Q$ with the same mesh that does not have $$|S(f, Q) - L| < 1.$$

Now $f$ must be unbounded on interval $[x_{j-1}, x_j]$ for some $j$. (If not, it's bounded on a finite collection of intervals, hence bounded.) Note that $j$ here is the index of just one interval. (I'll use "$k$" for a generic interval index.)

Therefore, for any $M$, there's a point $s \in [x_{j-1}, x_j]$ with $|f(s)| > M$.

I'll choose $$M = f(t_j) + \frac{10}{x_j - x_{j-1}}.$$ Let $s_j$ be the point with $|f(s_j)| > M$.

Now, let $Q$ be the partition with the same split-points as $P$, and the same tags as $P$, except replace the single tag $t_j$ with $s_j$. Let's compare $S(f, Q)$ to $S(f, P)$: $$S(f, P) - S(f, Q) = \sum_{k=1}^n f(t_k)(x_k−x_{k−1}) - \left[ \sum_{k=1}^{j-1} f(t_k)(x_k−x_{k−1}) + f(s_j) (x_j - x_{j-1}) \sum_{k = j+1}^n f(t_k)(x_k−x_{k−1}) \right]$$ All terms here cancel except for the $j$th, so \begin{align} |S(f, Q) - S(f, P) | &= | f(t_j)(x_k−x_{k−1}) - f(s_j) (x_j - x_{j-1}) |\\ &= |(f(t_j) - f(s_j))(x_j−x_{j−1}) | \\ &> |(\frac{10}{x_j - x_{j-1}})(x_j−x_{j−1}) | \\ &= 10 \end{align} Now since $|S(f, P) - L | < 1$ but $|S(f, Q) - S(f, P)| > 10$, the triangle inequality tells us that $|S(f, Q) - L | > 9$, and hence is not less than one.

• Thanks, this really helps Also I have already re-posted it following your suggestion so plesae ignore that as it will be a duplicate and I will delete it shortly. – user152874 Apr 5 '17 at 17:50

You never shared with us your definition of the Riemann integral. In Rudin, for example, boundedness of $f$ on $[a,b]$ is part of the definition.

I suppose you are letting $f:[a,b]\to \mathbb R,$ with no boundedness assumed, and are using the $\epsilon-\delta$ definition with Riemann sums. So $\int_a^b f =L$ means for every $\epsilon>0$ there is $\delta >0$ such that if $P$ is a partition with mesh size $<\delta,$ then

$$\tag 1 |\sum_{P}f(c_k)\Delta x_k -L| <\epsilon$$

for all choices of $c_k \in I_k.$ (Here the $I_k$ are the subintervals induced by $P,$ and $\Delta x_k$ is the length of $I_k.$)

Suppose $f$ is unbounded on $[a,b].$ WLOG $f$ is not bounded above. Take any partition $P$ of $[a,b].$ Then $f$ is not bounded above on some $I_{k_0}.$ This implies there is a sequence of points $d_m\in I_{k_0}$ such that $f(d_m) \to \infty.$ In the subintervals $I_k, k\ne k_0,$ take $c_k$ to be the left end point of $I_k.$ We'll keep the $c_k$'s fixed while letting $d_m$ vary. We get

$$\tag 2 \lim_{m\to \infty} \left (\sum_{k\ne k_0} f(c_k)\Delta x_k + f(d_m)\Delta x_{k_0}\right ) = \infty.$$

Now all Riemann sums in $(2)$ take place with respect to the partition $P.$ And this can be done for any $P.$ This destroys any hope that something like $(1)$ could hold.

format: plain


text:
# Evolutionary Games and Population Dynamics / Selection-diffusion / Gallery

Brave New World of Dynamical Pattern Formation Typical patterns for different efficiency of the public good $$r$$ and ratios of the diffusion coefficients. For $$r > r_\text{Hopf}$$ the co-existence equilibrium is stable and unstable otherwise, i.e. the population would invariably disappear in the absence of spatial extension. The spatial dynamics can be categorized into different kinds of pattern formation processes: homogenous distributions (yellow frame), diffusion induced instability - Turing patterns (green), diffusion induced co-existence (red), spatio-temporal chaos (blue) and extinction (black).
Effects of initial configurations Same as above but for symmetrical initial configurations with a disk in the center with equal densities of cooperators and defectors. The initial symmetry is often largely preserved by the deterministic dynamics but more importantly, the characteristic features of the emerging patterns remain unchanged.


text:
# Dynamics of matter-wave quantum emitters in a structured vacuum

@article{Stewart2020DynamicsOM,
title={Dynamics of matter-wave quantum emitters in a structured vacuum},
author={Michael Stewart and Joonhyuk Kwon and Alfonso Lanuza and D Schneble},
journal={arXiv: Quantum Gases},
year={2020}
}
• Published 5 March 2020
• Physics
• arXiv: Quantum Gases
The characteristics of spontaneous emission can be strongly modified by the mode structure of the vacuum. In waveguide quantum-electrodynamics based on photonic crystals, this modification is exploited to engineer atom-photon interactions near a band edge, but the physics of coupling to an entire band has not yet been explored in experiments. Here, we study the decay dynamics of an array of matter-wave quantum emitters coupled to a single band of an effective photonic crystal waveguide…
4 Citations
Undamped Rabi oscillations due to polaron-emitter hybrid states in a nonlinear photonic waveguide coupled to emitters
• Physics
Physical Review A
• 2022
The collective dynamics of two non-interacting two-level emitters, which are coupled to a structured wave guide that supports two-photon bound states, is investigated. Tuning the energy of the two
• Physics
• 2022
We propose a novel quantum electrodynamics (QED) platform where quantum emitters interact with a Hofstadter-ladder waveguide. We demonstrate several intriguing phenomena stemming from the nontrivial
Extensible quantum simulation architecture based on atom-photon bound states in an array of high-impedance resonators
• Physics
• 2021
Marco Scigliuzzo,1, ∗ Giuseppe Calajò,2 Francesco Ciccarello,3 Daniel Perez Lozano,1 Andreas Bengtsson,1 Pasquale Scarlino,4 Andreas Wallraff,5 Darrick Chang,2 Per Delsing,1 and Simone Gasparinetti1,

## References

SHOWING 1-10 OF 35 REFERENCES
Spontaneous emission of matter waves from a tunable open quantum system
• Physics
Nature
• 2018
This work demonstrates similar behaviour in a system of artificial emitters, realized using ultracold atoms in an optical lattice, which decay by emitting matter-wave, rather than optical, radiation into free space.
Markovian and non-Markovian dynamics of quantum emitters coupled to two-dimensional structured reservoirs
• Physics
• 2017
The interaction of quantum emitters with structured baths modifies both their individual and collective dynamics. In Gonz\'alez-Tudela \emph{et al} we show how exotic quantum dynamics emerge when QEs
Photon localization versus population trapping in a coupled-cavity array
• Physics
• 2014
We consider a coupled-cavity array (CCA), where one cavity interacts with a two-level atom under the rotating-wave approximation. We investigate the excitation transport dynamics across the array,
Photonic materials in circuit quantum electrodynamics
• Physics
• 2020
Photonic synthetic materials provide an opportunity to explore the role of microscopic quantum phenomena in determining macroscopic material properties. There are, however, fundamental obstacles to
Quantum matter built from nanoscopic lattices of atoms and photons
• Physics
Rochester Conference on Coherence and Quantum Optics (CQO-11)
• 2019
This Colloquium describes a new paradigm for creating strong quantum interactions of light and matter by way of single atoms and photons in nanoscopic lattices. Beyond the possibilities for
Quantum optics, what next?
• Physics
• 2017
Quantum optics is a well-established field that spans from fundamental physics to quantum information science. In the coming decade, areas including computation, communication and metrology are all
Spontaneous emission from a medium with a band spectrum
A theoretical analysis is made of the spontaneous emission from an excited atom or a molecule in a one-dimensional periodic structure formed by a screened two-conductor line filled with periodically
Cavity quantum electrodynamics
• Physics
• 2006
This paper reviews the work on cavity quantum electrodynamics of free atoms. In recent years, cavity experiments have also been conducted on a variety of solid-state systems resulting in many
Spontaneous and induced atomic decay in photonic band structures
• Physics
• 1994
Abstract We present a comprehensive quantum electrodynamical analysis of the interaction between a continuum with photonic band gaps (PBGs) or frequency cut-off and an excited two-level atom, which


text:
Mock modular forms and (indefinite) quadratic forms - MathOverflow most recent 30 from http://mathoverflow.net 2013-05-25T20:33:18Z http://mathoverflow.net/feeds/question/120324 http://www.creativecommons.org/licenses/by-nc/2.5/rdf http://mathoverflow.net/questions/120324/mock-modular-forms-and-indefinite-quadratic-forms Mock modular forms and (indefinite) quadratic forms Richard Eager 2013-01-30T14:12:03Z 2013-02-01T15:28:16Z
Define the function $$f(q,z,y) = \sum_{n \ge 0,m,l} c(n,m,l) q^n z^m y^l$$ where $c(n,m,l)$ is defined by $$c(n,m,l) = \begin{cases} (-1)^{s+l} & \mbox{if } 4n - m^2 + l^2 = 2s(s+1) \end{cases}$$ $$\begin{cases} 0 & otherwise \;\;\;\;\;\;\;\;\;\;\;\;\;\;
 \end{cases}$$ for some integer $s$ and $c(n,m,l) = 0$ unless $4n - m^2 -l^2 \ge 0.$ $f(q,z,1)$ is known to be related to a Mock modular form. I conjecture that $$f(q,1,-1) = \sum_{n \ge 0} (-1)^n (2n + 1) q^{n(n+1)/2}.$$ Is there an elementary proof of the above conjecture? Is the function $f(q,z,y)$ a known mathematical object, perhaps related to a Siegel modular form?

Update: $f(q,z,y)$ is a product of Jacobi theta functions and $\mu(q;z,y),$ where $\mu(q;z,y)$ is a Lerch sum studied by Zweger in his thesis. Zweger's thesis also relates mock modular forms to indefinite quadratic forms of signature $(1,n),$ so perhaps it isn't too unsurprising that $f(q,z,y)$ takes a "nice" form.


text:
# How to train a xgboost model on data that is too big for the memory?

What are the best practices to train xgboost (eXtreme gradient boosting) models on data that is to big to hold it in memory at once? Splitting the data and train multiple models? Are there more elegant solutions?

• Have you tried my suggestion below with the “boost from prediction”? Ironically, I now face similar problem and just wanted to check whether you were able to resolve it somehow. Thank you. Dec 21 '17 at 21:56
• Stacking would be a solution for sure. However, to make a prediction in the end you will have to run several models sequentially. It would be more elegant to generate subsets of a large dataset on the fly and build it into a single model. Jan 4 '18 at 8:12

You can train xgboost, calculate the output (margin) and then continue the training, see example in boost from prediction.

I‘ve not tried it myself, but maybe you could train on the first subset of your data (say 10%) and then continue on another subset, etc.

Update

Step by step procedure

1. Split the data into N manageable subsets, set n=1
2. Train xgboost on n-th subset
3. Calculate the prediction (margin) for n+1 subset using the model obtained from previous
4. Add the margin into the n+1 subset via setinfo
5. Increment n

Steps 2-5 to be repeated N times.

• i don't think this is what OP is asking. In your link it shows how to pass the WHOLE data through XGBoost and obtain new trees based using the scores from a model that was previous built. The OP is asking how to train the model on A PART of the data and then continue the training on another PART. Oct 24 '18 at 3:24

I don't think what you are asking for is possible. See this issue.

I understand that you want to train the model on A PART of the data and then continue the training on another PART and so. So @aivanov's answer will not help in this regard.

If you are using R, have you considered the bigmemory and ff packages?

I don't have much experience using these myself, but would be interested to see if they help with the issue at hand.


text:
# Direct summands of direct sum of line bundles on projective varieties

Serre-Swan's theorem (see the MO discussion) says that any locally free sheaf over an affine variety is a direct summand of a free sheaf. However, this is not true on projective varieties. It is not hard to check that a non-trivial line bundle with non-zero global sections can not be a direct summand of a free sheaf. The reason is that, being a direct summand of a free sheaf implies that the dual line bundle also has non-zero global sections. But that implies the line bundle is trivial.

I am wondering if the same holds for vector bundles, i.e. if a vector bundle and its dual over a projective variety both have non-zero global sections, then the vector bundle is trivial.

Another I think related question is the following:

Is a direct summand of a direct sum of line bundles on a projective variety also a direct sum of line bundles?

-

1) The vector bundle $E=\mathcal O(1) \oplus \mathcal O(-1)$ over $\mathbb P^1$ has non trivial sections and so has its dual (which happens to be the same bundle $E$ ).
However $E$ is non trivial because each of its global sections has a zero.

2) Yes, a direct summand of a direct sum of line bundles on a complete variety is also a direct sum of line bundles.
This follows from Atiyah's general version of the Krull-Schmidt theorem .

Edit
Since Fei asks, let me remark that the trick in 1) works also for bundles which are not direct sums of line bundles:
If $T$ is the tangent bundle to $\mathbb P^2$, then $E=T \oplus T^*$ has non-trivial global sections , and so has $E^*=E$
However, since $T$ is indecomposable $E$ is not a sum of line bundles by Atiyah's result and is thus a fortiori not trivial.

-
Thank you so much for your answer. For the first question, do you have a counterexample that $E$ is not a direct sum of line bundles? –  Fei YE Mar 25 '12 at 11:07
Dear @Fei,yes: I have added such a counterexample in an Edit. –  Georges Elencwajg Mar 25 '12 at 11:26
Dear George, Thank you so much for the example. I saw you used a fact $E^\ast=E$ for $E=T\oplus T^\ast$. Is there are general relation between $E$ and $E^\ast$. I saw somewhere that $E^\ast\otiems \det(E)=E$. Is this really true in general? Thanks again for you help! –  Fei YE Mar 25 '12 at 12:14
Dear @Fei, I have just used that the dual of a direct sum is the direct sum of the duals. The formula you quote is already false for every non trivial line bundle $L$ since $L^* \otimes det(L)=L^* \otimes L$ is trivial. –  Georges Elencwajg Mar 25 '12 at 12:24
Dear @Fei, yes, for rank $2$ it is correct: Hartshorne, Ex.II 5.16 –  Georges Elencwajg Mar 25 '12 at 14:06
show 1 more comment


text:
# Using .eps Graphics Files with PDF output

Oops — I forgot to thank the blog owner najmi for inviting me to be an author! How unbecoming of me! *grovel* Thanks again! 🙂

For quite some time now, the graphicx package would automatically sort out which file types to use whenever it encounters \includegraphics{some_file}, where the extension of some_file is not specified.

That is to say, if you are generating to a .dvi, the system would look for some_file.eps. If you are generating to a .pdf, it looks for some_file.{pdf|jpg|png}.

But what if you want a .pdf and you work with .pdf, .jpg and .png graphic files, but suddenly you received a .eps graphics from your collaborator? (This actually happened to a lecturer at my university.)

Well, you could always convert it manually to a .pdf yourself with GhostScript or GIMP, of course. But now the process can be automated if you’re using MikTeX 2.8 or TeXLive 2009, as highlighted by Rob Hyndman.

In MikTeX 2.8, use the epstopdf package together with graphicx, like so:

\documentclass{article}\usepackage{graphicx,epstopdf}\begin{document}\includegraphics[width=\textwidth]{some_file}\end{document}

If all you have is some_file.eps, it will be automatically converted to a some_file-eps-converted-to.pdf during the first time the LaTeX document is compiled and included in the final document output.

If you’re using TeXLive 2009, you don’t even have to load epstopdf explicitly, but I say include it anyway, for portability’s sake.

## 2 thoughts on “Using .eps Graphics Files with PDF output”

1. Great!

Just a bit wonder, why EPS is popularly used in LaTeX?

2. My guesses:

1. Historical reasons: In the beginning there were TeX and DVI and EPS. PDF came much, much later.

2. Dislike of proprietary formats (PDF).

3. PNGs and JPEGs are raster graphic formats and doesn't measure up to professional publishing requirements. EPS and PDF support vector graphics and therefore are preferred.

4. Generating and previewing DVIs are often faster than pdflatex and previewing PDFs (especially if you're using Acrobat's products).


text:
Running scientific instruments in Emacs and recording the results

| categories: | tags: | View Comments

Today we look at running a scientific instrument via http requests from Emacs and org-mode. We will use a Gamry Ref600 potentiostat because Gamry has very nicely provide a COM interface we can access via Python. This will be only a proof of concept to see what it is like. We will not consider any issues of security, etc…, only what is it like to do it.

The setup will look like this: we will run a flask web app that uses python to control the instrument via http requests. Why? Because I want to run the instrument from my Mac ;) and so far there are only Windows drivers for the instrument. So, we run the flask app on the Windows machine, and I run it from here on my Mac by sending requests. Flask takes care of converting requests to action using Python. You can see the Flask app here.

Let's see what is online:

curl jkitchin-win.cheme.cmu.edu:5000/pstats

(u'REF600-13089',)

We have one potentiostat online with serial number 13089. I have a dummy cell connected to it which has a little resistor on it. So we can run a cyclic voltammogram and it should be a straight line. We have to know a bit about what is returned. We will get a json file back, and it will have the data in it. The data will be a list of lists. The data we want is in columns 1 and 3 (python indexing). Obviously you need some prior knowledge of what data comes back to use this. That would come from reading some documentation.

import requests
import numpy as np
import matplotlib.pyplot as plt

resp = requests.get('http://jkitchin-win.cheme.cmu.edu:5000/cv?endv=0.25&startv=-0.25')

dj = resp.json()
data = np.array(dj['data'])

plt.plot(data[:, 1], data[:, 3])
plt.xlabel('Voltage (V)')
plt.ylabel('Current (A)')
plt.tight_layout()
plt.savefig('cv-1.png')

Well, there you have it. Possibly the first Gamry Ref600 to ever have been driven from a Mac ;) Let me be more explicit about that; I could also run this from Linux, an iPad, etc… You could do this in a browser, or in an IPython notebook, or in Matlab, among many other possibilities. You could write a script in perl, shell, ruby, emacs-lisp, or any other language that supports http requests.

I am not sure why the graph is not perfectly linear, maybe there is some capacitive charging that starts out. The resistance based on the current at 0.2V is about 2000 ohms, which is in good agreement with what is listed on the board the dummy cell is on.

1 Summary thoughts

There are a host of interesting issues one eventually has to consider here including security, but also error management and debugging. I hacked something like an http api here by running flask on the windows machine running the instrument. That is a layer of abstraction on an abstraction to start with. I think later instruments are likely to run these webservers themselves on small dedicated computers, e.g. via a Raspberry pi or Arduino chipset. It is not obvious how sophisticated you can make this with respect to triggering different instruments, etc…

In running this, my "notebook" was blocked while the experiment ran. It is possible to run things asynchronously, and sometimes that would make sense. In the example here, we have provided a very limited set of functions to "run" the potentiostat. It was only a proof of concept to get a sense for what it is like. In practice a fuller set of functions would be implemented. Another point to consider is how the data comes back from the potentiostat. We used json here because it is convenient, but we could just as well send files, and other sorts of data too.

This lays out the possibility to walk up to an instrument with an electronic notebook, setup and run the experiment, capture the results in the notebook and take it back to the office for analysis. Pretty cool.

So, here is my flask app. We setup a few routes using get requests to do things like get a list of the potentiostats online, and to run a cyclic voltamogram. As a side note, after this post is over, I am turning off the app, so you won't be able to repeat the codes ;) This is not a beautiful, secure or error tolerant code. It works enough for a proof of concept of simple experiments.

from flask import Flask, request, jsonify
import time

@app.route('/')
def hello_world():
return 'Hello World!'

@app.route('/pstats')
def get_pstats():
import win32com.client as client
devices = client.Dispatch('GamryCOM.GamryDeviceList')
result = str(devices.EnumSections())
return result

@app.route('/close_pstat')
def close():
import win32com.client as client
devicelist = client.Dispatch('GamryCOM.GamryDeviceList')

x = devicelist.EnumSections()[0]
pstat = client.Dispatch('GamryCOM.GamryPstat')
pstat.Init(x)

pstat.Close()

def run_ramp(Sinit,  # start value
Sfinal, # end value
ScanRate=1,
SampleRate=0.01,
CtrlMode=1,  # GamryCOM.PstatMode
fname=None):
'''We assume the first device is the one you want.
'''
import win32com.client as client
import numpy as np
devicelist = client.Dispatch('GamryCOM.GamryDeviceList')

x = devicelist.EnumSections()[0]

pstat = client.Dispatch('GamryCOM.GamryPstat')
pstat.Init(x)

pstat.Open()

dtaqcpiv=client.Dispatch('GamryCOM.GamryDtaqCpiv')
dtaqcpiv.Init(pstat)

sigramp=client.Dispatch('GamryCOM.GamrySignalRamp')
sigramp.Init(pstat, Sinit, Sfinal, ScanRate, SampleRate, CtrlMode)

pstat.SetSignal(sigramp)
pstat.SetCell(1) # 1 == GamryCOM.CellOn

try:
dtaqcpiv.Run(True)
except Exception as e:
pstat.Close()
raise

# NOTE:  The comtypes example in this same directory illustrates the use of com
# notification events.  The comtypes package is recommended as an alternative
# to win32com.
time.sleep(2) # just wait sufficiently long for the acquisition to complete.

acquired_points = []
count = 1
while count > 0:
count, points = dtaqcpiv.Cook(10)
# The columns exposed by GamryDtaq.Cook vary by dtaq and are
# documented in the Toolkit Reference Manual.
acquired_points.extend(zip(*points))

acquired_points = np.array(acquired_points)
if fname is not None:
np.savetxt(fname, acquired_points)

pstat.Close()
return jsonify({'data': acquired_points.tolist()})

@app.route('/cv')
def run_cv():
result = str(request.values)
startv = float(request.values.get('startv', -0.1))
endv = float(request.values.get('endv', 0.1))
scanrate = float(request.values.get('scanrate', 1.0))
samplerate = float(request.values.get('samplerate', 0.01))

data = run_ramp(startv, endv, scanrate, samplerate)
return data

if __name__ == '__main__':
app.run(host='jkitchin-win.cheme.cmu.edu', port=5000, debug=True)


text:
### Grant McDermott

Assistant Professor
Dept. of Economics
University of Oregon

# Water pricing in different sectors

My previous post linked to a new alarming/alarmist (potato, potahto) study of South Africa’s water pollution problems. I mentioned some economic elements that I feel need addressing if we’re ever going to provide the right set of incentives to ensure sustainable water use.

What I didn’t really discuss was the fact that water has a multitude of different uses, which can make valuation very hard. With regards to pricing structures, can you apply the same set of rules to all different users (residential vs agricultural vs industrial)? I can see more than a few practical problems, but there are certainly those who believe that we can and should do so. Prolific water blogger $$–$$ and recent addition to the job market! – David Zetland has long advocated a system of “all-in auctions”. In other words, very low prices for water aimed at meeting basic human needs, after which all comers are subject to compete against each other on price bids.

Having been raised in agricultural country, I can think of of a number of farmers that aren't going to be too happy with this type of proposal. However, that's a practical obstacle that can be overcome $$-$$ largely through greater dialogue and explanation of the actual procedures. Further to the point, a large number of farmers already control their own water sources in the form of small dams and so forth; you've obviously no need to submit to an auction if you already own the resource. In a similar vein, in instances where a bunch of farms are located along a river, the bidding could first be among that group of farmers. In cases where rivers flow to cities or industrial areas, up stream farmers could potentially make money by auctioning off "usage rights" to urban and industrial consumers.

On the subject of water usage by different consumers and sectors, I posted some relevant comments a while ago on a blog run by some Austrian friends (ideologically speaking that is, not geographically). With minor edits for style and following a discussion on the urgent need for a hike in (residential) water prices in SA:
One of the better articles on environmental issues published on this site. (Some of your climate and pollution arguments have been rather lacking in my opinion…) [Ed: For instance.]

Anyway, I certainly agree that water will be a defining issue for South Africa long before other crises and our current pricing structure is woefully inadequate. (As with many countries, the focus is on recuperating short-run average/operating costs instead of the economically sensible long-run marginal costs. While SRAC could be considered as constant, LRMC are typically increasing over time, so prices lie below their optimal value…)

However, the article does skirt around perhaps the most tricky issue with regards to water: There are MANY different uses for water and thus many different marginal values. Residential water consumption is fairly minute in comparison to agriculture and, to a lesser degree, industry. Thus $$–$$ and although it may not agree with your definition of inflation (i.e. an expansion of money supply) $$–$$ correct marginal pricing of water will certainly bring about an increase in general prices. Why? Water is a crucial factor of production for virtually everything we consume. Thus, an increase in price will carry through to all final products.

Not easy for society, but that is simply the bitter pill we need to swallow in order to establish sustainable use of what is perhaps our most valuable commodity.

After a fairly odd reply regarding my understanding inflation $$-$$ like all proud Austrians, they love to focus on inflation as a definitional good $$-$$ I ended an additional comment with:
[A]n increase in water prices will have much further implications than simply raised costs of showering or watering our gardens; it will affect prices for an entire host of economic goods.


text:
# diagonal var-cov matrix for random slope in lme4

I am fitting a mixed effect model with a random slope on a factor with 6 levels using the function lmer().

lmer(cc_marg ~ mPair + (0 + mPair | ratID), data = dat, REML = FALSE,
na.action = na.omit, control = lCtr)

I would like to simplify the variance-covariance matrix of the random effect by specifying that the 6 levels of the factor mPair are uncorrelated. In other words, I would like to specify a diagonal var-cov matrix. How can I do that using lmer()?

I know it would be easier to use lme() from the nlme package, as I would have access to pdClass. However, that package has issues at estimating the denominator DF for random slope models, and therefore, I would prefer to use the lme4 package. Thanks!