[WS2021] Homework 3 Discussion

[changkun]

This thread is open for the discussion of homework 3.

[torantie]

I have a problem :D

I am getting really high values for the linear equation and don't really know why (like -44104969.11230171). I am currently only using the Weight Matrix based on the code we received as a tip and omitted the mass matrix, which should be alright if i remember correctly. Up to let cholskyDecompositionMatrix = f.chol() i think everything is correct. What i am not sure about is the DenseMatrix let b = DenseMatrix.ones(numberOfVertices, 1) i use to solve the linear equation so maybe its the cause for the high values?

Current Code:

smooth(weightType, timeStep, smoothStep) {
      // TODO: implmeent the Laplacian smoothing algorithm.
      //
      try{
          // Hint:
          //
          //   1. Construct the Laplace matrix `L` for the given `weightType`

          let W = this.laplaceMatrix(weightType)
          let M = this.massMatrix(weightType)
          //   2. Solve linear equation: (I-tλL) f' = f using a Cholesky solver.
          let f = W.timesReal(timeStep)
          let cholskyDecompositionMatrix = f.chol()
          const numberOfVertices = this.vertices.length

          let b = DenseMatrix.ones(numberOfVertices, 1)
          let result = cholskyDecompositionMatrix.solvePositiveDefinite(b)



          //   3. Update the position of mesh vertices based on the solution f'.
          //
          for(let i = 0; i < numberOfVertices; i++){
              let resultEntryI = result.get(i,0)

              this.vertices[i].position = this.vertices[i].position.scale(resultEntryI)
          }
      } catch (e) {
          console.error(e)
      }
  }

/**
   * laplaceMatrix returns the Laplace matrix for a given laplaceType
   * @param {string} weightType indicates the type of the weight for
   * constructing the Laplace matrix. Possible value could be:
   * 'uniform', 'cotan'.
   */
  laplaceMatrix(weightType) {
    // TODO: construct the Laplace matrix.
    const numberOfVertices = this.vertices.length
    let weightTriplet = new Triplet(numberOfVertices, numberOfVertices)

    for (const vert of this.vertices) {
      const i = vert.idx
      let sum = 1e-8 // Tikhonov regularization to get strict positive definite

      vert.forEachHalfEdge(h => {
        let w = 0

        switch (weightType) {
          case 'uniform':
            w = 1
            break
          case 'cotan':
            w = (h.cotan() + h.twin.cotan())/2
        }

        sum += w

        weightTriplet.addEntry(w, i, h.twin.vertex.idx)
      })
      weightTriplet.addEntry(-sum, i, i)
    }

    let weightMatrix = SparseMatrix.fromTriplet(weightTriplet)

    return weightMatrix
  }

  massMatrix(weightType){
    const numberOfVertices = this.vertices.length
    let massTriplet = new Triplet(numberOfVertices, numberOfVertices)

    for (const vert of this.vertices) {
      const i = vert.idx

      switch (weightType) {
        case 'uniform':
          let neighbours = 0
          vert.forEachHalfEdge(() => {
            neighbours++
          })

          massTriplet.addEntry(neighbours,i,i)
          break
        case 'cotan':
          massTriplet.addEntry(vert.calculateVoronoiArea(),i,i)
      }
    }

    let massMatrix = SparseMatrix.fromTriplet(massTriplet)

    return massMatrix
  }

[changkun]

let b = DenseMatrix.ones(numberOfVertices, 1)

This is definitely a mismatch of the Laplacian smooth. You need to think more about the backward Euler that we use for time descritization. Under the Laplacian smooth case, the right hand side is the position of the vertices, therefore you will have to construct the b similar as follows:

let b = DenseMatrix.zeros(this.vertices.length, 3)
for (const v of this.vertices) {
   b.set(v.position.x, v.idx, 0)
   b.set(v.position.y, v.idx, 1)
   b.set(v.position.z, v.idx, 2)
}
const x = cholskyDecompositionMatrix.chol().solvePositiveDefinite(b)

Remember: If you get anything weird, it is most likely because of "garbage in garbage out", if you use the cotan form, you must have to verify it is a valid solution (recall the properties of cotan Laplace and when you can use Cholesky solver), i.e. the solution x that Cholesky solver gave you indeed satisfy Ax=b.

[torantie]

Ah okay i wasn't really sure what "b" was supposed to be. Makes sense that we need three coordinates to shift the vertices.^^
Unfortunately the values of the linear equation are still extremely high and deform the mesh. does this have something to do with the smooth step? I couldn't really understand how the step in the Laplacian Smoothing is used in the Algorithm.
But since the Default of the Smooth step is 1 i don't think it's the reason for these high numbers.
Is it perhaps wrong to omit the Mass Matrix?
Sorry for all the questions i am just really stuck^^

smooth(weightType, timeStep, smoothStep) {
    // TODO: implmeent the Laplacian smoothing algorithm.
    //
    try{
      // Hint:
      //
      //   1. Construct the Laplace matrix `L` for the given `weightType`

      let W = this.laplaceMatrix(weightType)
      let M = this.massMatrix(weightType)
      //   2. Solve linear equation: (I-tλL) f' = f using a Cholesky solver.
      let f = W.timesReal(timeStep)
      let cholskyDecompositionMatrix = f.chol()
      const numberOfVertices = this.vertices.length

      let b = DenseMatrix.zeros(numberOfVertices, 3)

      for (const v of this.vertices) {
        b.set(v.position.x, v.idx, 0)
        b.set(v.position.y, v.idx, 1)
        b.set(v.position.z, v.idx, 2)
      }
      let result = cholskyDecompositionMatrix.solvePositiveDefinite(b)

      //   3. Update the position of mesh vertices based on the solution f'.
      //
      for(let i = 0; i < numberOfVertices; i++){
        let resultX = result.get(i,0)
        let resultY = result.get(i,1)
        let resultZ = result.get(i,2)
        let currentPosition = this.vertices[i].position
        currentPosition.x = currentPosition.x * (resultX)
        currentPosition.y = currentPosition.y * (resultY)
        currentPosition.z = currentPosition.z * (resultZ)
      }
    } catch (e) {
      console.error(e)
    }
  }

[changkun]

Your implementation towards the right direction, but needs a little bit help with the mesh structure.

You might need to notice: this.vertices[i].position is permanently changed after one step of smooth because you updated them in the smooth function. This is what I mentioned in the class where you need a little bit change to the mesh structure to cache the original vertex position, and always use them to solve the Laplace equation, then update this.vertices for rendering the smoothed positions. Reading main.js might be helpful to understand the reason.

[torantie]

ah so the smooth step just defines how often the smooth algorithm has to be executed. Meaning i need something like this?:

smooth(weightType, timeStep, smoothStep) {
for (let j = 0; j< smoothStep; j++) {
      // TODO: implmeent the Laplacian smoothing algorithm.
      //
      try {
        // Hint:
        //
        //   1. Construct the Laplace matrix `L` for the given `weightType`
        let W = this.laplaceMatrix(weightType)
        let M = this.massMatrix(weightType)
        //   2. Solve linear equation: (I-tλL) f' = f using a Cholesky solver.
        let f = W.timesReal(timeStep)
        let cholskyDecompositionMatrix = f.chol()
        const numberOfVertices = this.vertsOrig.length
        let b = DenseMatrix.zeros(numberOfVertices, 3)

        for (const v of this.vertsOrig) {
          b.set(v.position.x, v.idx, 0)
          b.set(v.position.y, v.idx, 1)
          b.set(v.position.z, v.idx, 2)
        }

        let result = cholskyDecompositionMatrix.solvePositiveDefinite(b)

        //   3. Update the position of mesh vertices based on the solution f'.
        //
        for (let i = 0; i < numberOfVertices; i++) {
          let resultX = result.get(i, 0)
          let resultY = result.get(i, 1)
          let resultZ = result.get(i, 2)
          let currentPosition = this.vertices[i].position
          currentPosition.x = currentPosition.x * (resultX)
          currentPosition.y = currentPosition.y * (resultY)
          currentPosition.z = currentPosition.z * (resultZ)
        }
      } catch (e) {
        console.error(e)
      }
    }
  }

[changkun]

Bingo :)

[changkun]

Close due to the end of the submission period.