• This topic has 53 replies, 9 voices, and was last updated 4 years ago by YY.
Viewing 9 reply threads
  • Author
    Posts
    • #6168 Score: 0
      Scott
      Participant
      -3 pt

      目前完成情况是显示欧拉不会收敛,需要用上课时说的改进欧拉让他可以收敛吗?(如果减小delta t,飞出去需要的时间久一点),半隐式欧拉比较正常,主要问题是欧拉方式里的整体阻尼是如何实现,如果用1.4里的那个阻尼公式是有比较好的效果,但是这应该是用在Verlet方法里的所以感觉不对。
      另外pdf上描述的Verlet不是很理解,烦请大佬帮忙解释一下。

    • #6169 Score: 0
      Bear
      Participant
      -11 pt

      我和你在纠结一样的问题。。。
      – 欧拉法的全局 dumping 我把上课讲得两种方法都实现了,但最终会变成有惯性的绳子摆来摆
      去还是停不下来。 用 Verlet 的方法倒是能停
      – Verlet 看 wiki 就比较清楚,实际上就是通过一阶泰勒推出了下一时间点的位置可以通过当前时间点位置和上一时间点位置进行计算
      – 问题在于初始时怎么算上一时间点(还没走了),
      https://wenku.baidu.com/view/48619b136edb6f1aff001f51.html
      讲了两种方法,结果还是很奇怪直接掉下去了。。。
      – 而且还照上述链接的方法计算的话我就不知道在 Spring 里要算啥了。。。约束的实现不
      是在 mass 的 loop 里实现了?

      • #6171 Score: 0
        Scott
        Participant
        -3 pt

        课上讲的Adaptive Step Size用在整体阻尼上吗

        • #6172 Score: 0
          Bear
          Participant
          -11 pt

          ??? 那些不都是显式欧拉法的改进,和阻尼没关系?

          和阻尼有关系的在倒数第二节

      • #6178 Score: 1
        JZZ
        Keymaster
        2 pts

        我在半隐式欧拉的Add global damping位置加了一句:
        m->velocity = (1 – 0.00005)*m->velocity;
        就可以停下来了。因为感觉阻尼的存在就是和速度对抗,让速度降下来(不过老师讲要根据弹簧两个端点的相对速度计算,这么做感觉有点偷懒了。。。真是迷)
        显示欧拉,因为一开始速度是NULL,就在Mass的构造函数给velocity设了初值(0,0)(顺手把last_position也给了(0,0)的初值)然后使用这一帧的速度计算当前位置,结果绳子直接起飞了。。。。。真是不知道怎么写才对了。。。。
        verlet可以用当前质点位置的前一帧、后一帧的二阶泰勒展开解出来,不过里面带了一个加速度项;那求加速度不是又得在物理上分析受力了吗?感觉没有啥思路实现。。。。
        感觉这次实验写的好迷啊,哪都不清楚,求点拨

        This post has received 1 vote up.
        • #6183 Score: 0
          Mxkkk
          Participant
          1 pt

          那个加速度项直接用重力加速度就好,因为解约束的话就不考虑弹力了

          • #6186 Score: 0
            o_o_o_o_o
            Participant
            -1 pt

            插楼同问这里的 global damping kd取值,我的理解是跟你差不多,这里global damping 就是取一个 f = -kd b_dot

            • #6187 Score: 0
              Bear
              Participant
              -11 pt

              他直接改了速度了。。。公式不是应该算的是力?

              • #6195 Score: 0
                JZZ
                Keymaster
                2 pts

                计算阻尼力的话,感觉在springs循环合适吧。因为要用弹簧两端的质点计算相对速度,但题目提示在masses遍历里加一个global的阻尼。。。就不知道怎么理解了。话说大大Vertel方法怎么实现呐,加速度取gravity的话绳子就从一端落下去了,或者是我springs里的长度约束没写对?

                • This reply was modified 4 years ago by JZZ.
                • #6197 Score: 0
                  Bear
                  Participant
                  -11 pt

                  – 阻尼有两种方法啊,一种在 Spring 一种在 Mass。我在下面都写了,不过这两种我发现都没用。。。
                  – Vertel 的话你忘了固定的 mass 是不能动的?Spring 里要注意固定的不要动

                  • #6198 Score: 0
                    JZZ
                    Keymaster
                    2 pts

                    学到了,学到了。不过阻尼的话我也在springs里计算过,现在想想在springs里计算然后算到合力里,其实就相当于变相改变弹力了;在masses里计算阻尼再去改变合力,跑回来一圈感觉又是变相的改变重力,会继续参加到下一步的速度计算(整个系统感觉没有损耗,所以能一直动)
                    如果直接对速度打折的话,相当于阻尼作用结果的直接使用,就能从系统中不断损耗能量

                    • This reply was modified 4 years ago by JZZ.
                    • #6200 Score: 0
                      Bear
                      Participant
                      -11 pt

                      ???
                      一个物体运动受到多个方向不同的力没问题啊,每次给的阻尼也都是反向的力会减少整体系统中的能量。

                      我觉得你这个例子不对劲,而且如果单纯打对折有用完全正确的话显式欧拉的方法最终就不会飞了。。。

                      • #6203 Score: 0
                        JZZ
                        Keymaster
                        2 pts

                        好久没学物理,脑袋废了;但可以把节点设置成两个跑一下,应该能看出直接对力修正是没有办法影响摆动高度的,应该是这个框架不是严格的物理仿真?或者是我写法有问题。。。

                      • #6204 Score: 0
                        Bear
                        Participant
                        -11 pt

                        反正在我这里是两种方法都停不下来;其中在 Spring 计算的方法是能变成正常幅度的钟摆然后停不下来;你可以看下我下面的代码有问题一起讨论讨论

        • #6184 Score: 0
          Bear
          Participant
          -11 pt

          – 全局阻尼那个,TODO 在 Mass 处而不是 Spring 处,所以可能是按最普通的来不计算内力。我是奇怪我用两端点算最后结果是根据惯性定律摆来摆去停不下来。。。
          m->forces += -0.00005 * m->velocity.unit(); 最普通的这样也不行
          – Vec 初值就是 0 啊,last_pos 初始值和 pos 一样
          – 加速度是因为有重力啊,我是看 Verlet 本身的力有说根据势能算,不过按 PDF 的说法应该就是按我上面说得在 Spring 强行移动恢复原长, Mass 再更新位置

          • #6189 Score: 0
            Bear
            Participant
            -11 pt

            我觉得我实现也没问题啊🤔

            Attachments:
            You must be logged in to view attached files.
            • #6191 Score: 0
              o_o_o_o_o
              Participant
              -1 pt

              你的是没问题,我觉得这算 Internal Damping for Spring, cpp文件里 hint 不是要求添加 global damping,我理解的global damping 是针对每个质点 -kd * 物体速度

              • #6192 Score: 0
                Bear
                Participant
                -11 pt

                – 对就是这个代码最终 spring 也停不下来,变成了惯性的正常幅度运动。。。
                – 我也是这么想的,那直接不算内力算 damping 应该是

                
                  m->forces += -0.00005 * m->velocity.unit();
                

                直接算速度的反方向加到力上,但实际上也没用我发现。。。

                • #6193 Score: 0
                  o_o_o_o_o
                  Participant
                  -1 pt

                  对, 你是对的,这样加 global damping 是停不下来的, o(╯□╰)o

                  • #6194 Score: 0
                    Bear
                    Participant
                    -11 pt

                    所以作业就很迷。。。两种做法都没用

    • #6173 Score: -1
      Bear
      Participant
      -11 pt

      Verlet 在 Mass 处要移动点使得弹簧维持原长,虽然我这里实现起来好像变成了直挺挺的。。。

      This post has received 1 vote down.
      • #6174 Score: 0
        Scott
        Participant
        -3 pt

        那他这个Spring处要干啥呀?

        • #6175 Score: 0
          Bear
          Participant
          -11 pt

          就是像我说的改位置,在 Mass 处再根据给的公式移动到下一个位置(大概)

    • #6176 Score: 0
      安安
      Participant
      -1 pt

      请问一下Verlet的x(t-1)是怎么得到的呢

    • #6177 Score: 0

      我理解的x(t-1)是那个m->last_position,也不知道是不是的

      • #6185 Score: 0
        Bear
        Participant
        -11 pt

        对。然后记得维护更新,

        初始 0 时的 t-1 可能要靠近似?

        • #6188 Score: 0
          Bear
          Participant
          -11 pt

          试了下初始 0 时的 t-1 不用写近似也行

    • #6206 Score: 0
      Scott
      Participant
      -3 pt

      所以说对于显式欧拉存在方法让它不飞出去吗,以及对于半隐式欧拉在mass处的global damping除了单纯地降低速度没有什么有效的方法可以让它停下来吗?

      • #6207 Score: -1
        Bear
        Participant
        -11 pt

        。。。你还是再看看课吧,概念都没搞懂

        This post has received 1 vote down.
        • #6208 Score: 0
          Scott
          Participant
          -3 pt

          我看课只讲了弹簧内力,在mass就只能用最简单粗暴的方法?

          • #6209 Score: 0
            Bear
            Participant
            -11 pt

            内力和简单的那种我都有写发现没啥效果。。。
            我代码在上面你可以看看我是不是打错了

    • #6210 Score: 0
      Atlas
      Participant

      http://datagenetics.com/blog/july22018/index.html
      我看了这个网站写的verlet,看起来结果还蛮正常的。
      1. 遍历所有mass,加速度用重力,算出新的position。
      2. 遍历所有string,根据第一步的两边mass的position算出实际长度,然后要让它恢复到rest_length,以此更新两边mass的position(注意是否pinned)。

      • This reply was modified 4 years ago by Atlas.
      • #6212 Score: 0
        Bear
        Participant
        -11 pt

        这不是一码事吗。。。

        框架最开始时Spring 长度没有发生变化,所以第一步时的 Spring 中实际上啥也没干

    • #6213 Score: -1
      Atlas
      Participant

      抱歉,前面回复太多了,我可能没看懂dalao写的;;
      Verlet那里,spring和mass都有用。spring里在根据弹簧长度constraint修正点的位置。mass里在用公式算点的位置。
      Euler阻尼的话,我是在mass那里加m->forces += -0.005 * m->velocity;,就能停下来。

      This post has received 1 vote down.
      • #6214 Score: 0
        Bear
        Participant
        -11 pt

        欧拉法这里显式隐式用这个都能停?
        那可太神奇了。。。

      • #6215 Score: 0
        Bear
        Participant
        -11 pt

        请问我这些写的和你差别很大?
        我感觉我这里不应该有问题啊😂

        Attachments:
        You must be logged in to view attached files.
    • #6218 Score: 0
      Atlas
      Participant

      忘记说了,显式的停不下来。
      我是把m->forces += -0.005 * m->velocity;加在算加速度之前的,你这样加在最后,不是就没用嘛。因为后面还有一句m->forces = Vector2D(0, 0);

      • #6221 Score: -1
        Bear
        Participant
        -11 pt

        有道理。。。我被TODO影响了

        不过在 Spring Loop 里算内力的那种也停不下来啊。。。

        This post has received 1 vote down.
        Attachments:
        You must be logged in to view attached files.
        • #6249 Score: 0
          YY
          Participant
          1 pt

          阻尼与速度乘反比f=-kd*x’,你把Kd放在spring中,实际上是在削减弹簧的刚度ks,刚度不断衰减最多就是刚度为0,最终只受重力,自由落体;此外,也很难保证在速度方向上的衰减,加速度是弹力和重力的合力,速度是由加速度得到的,这样质点速度大概率与弹簧并不在同一条直线上,更不用说方向相反了。还有一个疑问,框架里的gravity到底是重力加速度还是重力啊,有点分不清了,翻译是重力,但是看见有人乘mass。

          • #6251 Score: 0
            Bear
            Participant
            -11 pt

            我在spring 里算是因为用考虑内力的方法

            如果在mass里没法算内力啊

            应该是重力加速度,因为是config定义的,我当成g也做完了

            • #6256 Score: 0
              YY
              Participant
              1 pt

              我输出了一下,mass和gravity大小都是1,好像当成什么都一样。
              力的作用是加速度,然后影响速度,最后反应到位移上。所以不论是先算合力,然后算速度位移,还是先分别计算每个分力的速度位移再相加最后都是一样的。
              所以可以在TODO dumping中用v*kd计算阻力→加速度→速度增量→速度。因为欧拉法本来就是近似,而且每个步长速度增量很小,所以放在哪里计算都一样;不过在masses中计算更简单,因为阻力实际是作用在质点上的,如果放在springs中一个弹簧两个质点,相邻弹簧由公共质点,可能计算重复。
              如果按照老师讲的阻尼取决于相对速度的话,在springs中计算更方便,但是我认为这个地方取绝对速度也可以,我以前上结构动力学的时候,粘滞阻尼是与运动方向相反,这里的运动方向都是基于地面(地球),而不是质点相对速度;因为阻尼最终导致运动停止,基于绝对速度可以达到这个效果(我记得阻尼系数是根据实验得来的,这里的实验应该是绝对速度,实验应该不会还分出各个质点求相对速度)。基于相对速度,各部分相对静止,最后也是基于地面静止,效果相同。我的理解是基于绝对速度,阻尼可以当作空气阻力,基于相对速度可以看出能量以热量耗散来损失。
              还有那个显示欧拉,我把Application::render() 里调用simulateEuler循环改成了i<1,速度放慢了,看见不是直接发散,而是正常运动了一会发散了,应该是老师讲的稳定性问题。

              • #6257 Score: 0
                Bear
                Participant
                -11 pt

                辛苦了这么多字。。。

                – F = -kd*V 算阻力自然还是放在 Mass Loop 好
                – 我指得是 F = -kd * dot(difference(v2, v1), difference(pos2, pos1).unit()) 这个利用内力算阻力的放在 Spring Loop 里才能算,同时这种写法没起到停止的作用(助教发帖这个不用实现了,估计可能框架没支持这种方法算)
                – 显式欧拉我靠把每帧 step 加到 128000 才最终稳定停下来,你可以实验试试

                • #6273 Score: 0
                  YY
                  Participant
                  1 pt

                  我之前把你在springs里的代码看错了。。。我用你这种阻尼加半隐式方程,发现最后会发散,你的f_d_ab和f_d_ba好像犯了。我反过来之后不发散单也不会使得运动停止。我是这么理解的:
                  我把节点数定位2,一端固定就成了钟摆运动。质点受到重力,拉力,阻力,其中阻力这项按照ppt计算:大小取决于相对速度,方向为相对速度沿弹簧方向的投影。实际只有两个方向的力:一个为弹簧方向的合力(阻力,拉力,重力分量)改变速度方向,一个为速度方向(重力分量),这样一来速度方向的加速度基本始终不变(每个周期不变)。假设由于弹簧伸缩,不是精确的圆周运动,将运动分解,一部分为圆周运动,一部分为径向运动,这样我们的阻力始终只影响径向运动。从能量守恒的角度,阻力为绳子径向,运动方向几乎垂直于径向,径向位移很小,阻力也很小,合起来的做工就更能小了,而且按照这种算法,这部分能量不一定是真的损失了,因为这个值可能正负交替。
                  我回去又查了查这个能量耗散,一种就是类似空气阻力的外部分损失cv’,一种是内部能量耗散发热,这部分跟材料有关,主要取决于应力应变曲线(可以当作位移和力的曲线),斜率就是刚度,所围面积为做功。如果加载和卸载都是沿同一条曲线,那么就没有能量耗散,实际并不是这样,举其中一个情况:发塑形变形,当卸力时,只发生弹性恢复,塑形形变不回复,这部分就是耗散的能量。
                  所以我觉得,老师上课讲的那种基于径向相对速度得到的阻力本身可能就不正确,从改变弹簧刚度的角度讲,无法区分卸载和加载。从速度的方向讲在质点运动方向做功微乎其微。

                  Attachments:
                  You must be logged in to view attached files.
                  • #6275 Score: 0
                    YY
                    Participant
                    1 pt

                    上升沿曲线面积就是加载向弹簧做的功,下降沿曲线就是卸载时弹簧向外做的功,差值为能量耗散。

                    • #6276 Score: 0
                      Bear
                      Participant
                      -11 pt

                      – 图形学的模拟和物理仿真目标又不一样,公式无法反应真实世界很正常。Phong 那么假但现在不还用吗
                      – 但这个系统已经诞生很久了,早年间的布料模拟一直用这个,至少说明能停下来
                      – 至于现在算内力停不下来应该是有其它因素框架并没有考虑
                      – 我的打法有没问题也没法验证,但按照公式应该是没错的。反正我的结果是变成了正常的钟摆摆来摆去,幅度较小

                      • #6278 Score: 0
                        YY
                        Participant
                        1 pt

                        我的意思是,如果按照课件公式的话你的f_d_ab不应该是fb吗,或者我把课件公式理解错了。
                        还有就是我很认同图形学不是仿真,所以在哪里damping,用什么方法算阻尼能够正常的停下来就行了。
                        我以上只是对为什么停不下来的猜测,没有否定图形学的研究成果的意思。布料能停下来应为每个节点处在网络中而不是线中,多了两个约束,这个方法让相邻节点径向保持相对静止,但是可以圆周运动,不相邻的节点就有相对运动。

        • #6250 Score: 0
          YY
          Participant
          1 pt

          矢量具有可加性,其实不论在springs中还是在masses中算都是一样的,只不过在springs中算的是力,这个力要用-kd*v得到;在masses中算的是力作用的结果,而不能算力,因为最后力清零了。

          • #6252 Score: 0
            Bear
            Participant
            -11 pt

            所以框架把TODO dumpling 放在最后面是迷惑行为?

    • #6226 Score: 0
      Atlas
      Participant

      我按照你那种方法算隐式内力,kd设成0.1就能停下来。kd更小的时候,最后就会看起来在oscillate。感觉已经是某些控制系统的东西了┑( ̄Д  ̄)┍

      • #6228 Score: 0
        Bear
        Participant
        -11 pt

        我设成0.1显式和隐式都停不下来 ╮(╯▽╰)╭
        直接就飞了。。。

        明明公式很简单我也没打错
        所以太难了。。。

        • This reply was modified 4 years ago by Bear.
        • #6231 Score: 0
          Atlas
          Participant

          是我写错了。。。我以为*和dot一样,结果*被重载成complex的计算了。
          所以,显式的怎样都停不下来,隐式的只能在mass那里加kd*v才能停下来。
          我不说了,还是等dalao们了,要是能删帖就好了;;

          • #6232 Score: 0
            Bear
            Participant
            -11 pt

            所以我已经怀疑人生了。。。

            看看后来人能不能躺平路吧,唉

          • #6240 Score: 1
            Bear
            Participant
            -11 pt

            还是因为显式欧拉不稳定, 我把每帧步骤调到了 12800,终于显式欧拉也停下来了

            This post has received 1 vote up.
            • #6259 Score: 0
              o_o_o_o_o
              Participant
              -1 pt

              隔壁有个hw8的解答贴,发现如果只加 global damping 的话显示欧拉是会停下来的,kd的大小控制‘收敛’的速度,你把kd调大一点,比如0.1看看,说不定就觉得好看了

Viewing 9 reply threads
  • You must be logged in to reply to this topic.