Mercurial > repos > pointsource_algs / file comparison
comparison: src/kernels/base.rs
src/kernels/base.rs
- changeset 72
- e9a460a0e638
- parent 61
- 4f468d35fa29
equal
deleted
inserted
replaced
| 71:e2953ffd4e0b | 72:e9a460a0e638 |
|---|---|
| 435 lip: F, | 435 lip: F, |
| 436 dbound: F, | 436 dbound: F, |
| 437 dlip: F, | 437 dlip: F, |
| 438 ) -> DynResult<F> { | 438 ) -> DynResult<F> { |
| 439 // For arbitrary ψ(x) = ∏_{i=1}^n ψ_i(x_i), we have | 439 // For arbitrary ψ(x) = ∏_{i=1}^n ψ_i(x_i), we have |
| 440 // ψ(x) - ψ(y) = ∑_i [ψ_i(x_i)-ψ_i(y_i)] ∏_{j ≠ i} ψ_j(x_j) | 440 // ψ(x) - ψ(y) = ∑_i [ψ_i(x_i)-ψ_i(y_i)] ∏_{j < i} ψ_j(y_j) ∏_{j > i} ψ_j(x_j) |
| 441 // by a simple recursive argument. In particular, if ψ_i=g for all i, j, we have | 441 // by a simple recursive argument. In particular, if ψ_i=g for all i, j, we have |
| 442 // |ψ(x) - ψ(y)| ≤ ∑_i L_g M_g^{n-1}|x-y|, where L_g is the Lipschitz factor of g, and | 442 // |ψ(x) - ψ(y)| ≤ ∑_i L_g|x_i-y_i| M_g^{n-1} ≤ √n L_g M_g^{n-1} |x-y|₂ |
| 443 // M_g a bound on it. | 443 // where L_g is the Lipschitz factor of g, and M_g a bound on it. |
| 444 // | 444 // |
| 445 // We also have in the general case ∇ψ(x) = ∑_i ∇ψ_i(x_i) ∏_{j ≠ i} ψ_j(x_j), whence | 445 // We also have in the general case [∇ψ(x)]_i = ψ_i'(x_i) ∏_{j ≠ i} ψ_j(x_j), whence |
| 446 // using the previous formula for each i with f_i=∇ψ_i and f_j=ψ_j for j ≠ i, we get | 446 // from above, the Lipschitz factor of [∇ψ(x)]_i is |
| 447 // ∇ψ(x) - ∇ψ(y) = ∑_i[ ∇ψ_i(x_i)∏_{j ≠ i} ψ_j(x_j) - ∇ψ_i(y_i)∏_{j ≠ i} ψ_j(y_j)] | 447 // L' = √(L_{g'}^2M_g^{2(n-1)} + (n-1)L_g^2M_{g'}^2M_g^{2(n-2)}) if n > 2, |
| 448 // = ∑_i[ [∇ψ_i(x_i) - ∇ψ_j(x_j)] ∏_{j ≠ i}ψ_j(x_j) | 448 // in particular |
| 449 // + [∑_{k ≠ i} [ψ_k(x_k) - ∇ψ_k(x_k)] ∏_{j ≠ i, k}ψ_j(x_j)]∇ψ_i(x_i)]. | 449 // L' = √(L_{g'}^2M_g^{2(n-1)} + L_g^2M_{g'}^2) if n > 2, and |
| 450 // With $ψ_i=g for all i, j, it follows that | 450 // Now the Lipschitz factor of ∇ψ is √n L' from |
| 451 // |∇ψ(x) - ∇ψ(y)| ≤ ∑_i L_{∇g} M_g^{n-1} + ∑_{k ≠ i} L_g M_g^{n-2} M_{∇g} | 451 // ‖∇ψ(x)-∇ψ(y)‖ = √(∑_{i=1}^n [∇ψ(x)-∇ψ(y)]_i^2) ≤ √(∑_{i=1}^n (L')^2) ≤ √n L' |
| 452 // = n [L_{∇g} M_g^{n-1} + (n-1) L_g M_g^{n-2} M_{∇g}]. | 452 if N > 2 { |
| 453 // = n M_g^{n-2}[L_{∇g} M_g + (n-1) L_g M_{∇g}]. | |
| 454 if N >= 2 { | |
| 455 Ok(F::cast_from(N) | 453 Ok(F::cast_from(N) |
| 456 * bound.powi((N - 2) as i32) | 454 * (dlip.powi(2) * bound.powi(2 * (N - 1) as i32) |
| 457 * (dlip * bound + F::cast_from(N - 1) * lip * dbound)) | 455 + F::cast_from(N - 1) |
| 456 * lip.powi(2) | |
| 457 * dbound.powi(2) | |
| 458 * bound.powi(2 * (N - 2) as i32)) | |
| 459 .sqrt()) | |
| 460 } else if N == 2 { | |
| 461 Ok((F::TWO * ((dlip * bound).powi(2) + (lip * dbound).powi(2))).sqrt()) | |
| 458 } else if N == 1 { | 462 } else if N == 1 { |
| 459 Ok(dlip) | 463 Ok(dlip) |
| 460 } else { | 464 } else { |
| 461 Err(anyhow!("Invalid dimension")) | 465 Err(anyhow!("Invalid dimension")) |
| 462 } | 466 } |
